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Abstract 

The  ability  to  locate  an  RF  transmitter  is  a  topic  of  growing  interest  for  civilian  and 
military  users  alike.  Geolocation  can  provide  critical  information  for  the  intelligence 
community,  search  and  rescue  operators,  and  the  warfighter.  The  technology  required  for 
geolocation  has  steadily  improved  over  the  past  several  decades,  allowing  better 
performance  at  longer  baseline  distances  between  transmitter  and  receiver.  The  expansion 
of  geolocation  missions  from  aircraft  to  spacecraft  has  necessitated  research  into  how 
emerging  geolocation  methods  perform  as  baseline  distances  are  increased  beyond  what 
was  previously  considered.  The  Cubesat  architecture  is  a  relatively  new  satellite  form 
which  could  enable  small-scale,  low-cost  solutions  to  USAF  geolocation  needs. 

This  research  proposes  to  use  Cubesats  as  a  vehicle  to  perform  geolocation  missions 
in  the  space  domain.  The  Cubesat  form  factor  considered  is  a  6-unit  architecture  that 
allows  for  6000  cm^  of  space  for  hardware.  There  are  a  number  of  methods  which  have 
been  developed  for  geolocation  applications.  This  research  compares  four  methods  with 
various  sensor  configurations  and  signal  properties.  The  four  methods’  performance  are 
assessed  by  simulating  and  modeling  the  environment,  signals,  and  geolocation 
algorithms  using  Matlab. 

The  simulations  created  and  run  in  this  research  show  that  the  angle  of  arrival  method 
outperforms  the  instantaneous  received  frequency  method,  especially  at  higher  SNR 
values.  These  two  methods  are  possible  for  single  and  dual  satellite  architectures.  When 
three  or  more  satellites  are  available,  the  direct  position  determination  method 
outperforms  the  three  other  considered  methods. 
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RADIO  FREQUENCY  EMITTER 
GEOEOCATION  USING  CUBESATS 


I.  Introduction 

The  problem  of  RE  transmitter  geolocation  has  been  a  topic  of  interest  for  military  and 
commercial  applications  since  shortly  after  World  War  IT  The  technology  began 
with  purely  land-based  apparatus,  designed  to  find  the  direction  of  arrival  by  comparing 
the  power  of  received  signals  among  apertures  in  antenna  arrays.  The  technology  has 
since  matured  and  been  applied  in  more  diverse  areas,  especially  in  the  last  several 
decades  as  communication  technology  has  grown  more  rapidly.  The  geolocation  platforms 
have  also  expanded  from  land-based  to  sea-  and  air-based,  and  are  beginning  to  be  applied 
in  space-based  applications. 

Applications  for  geolocation  techniques  are  wide  ranging,  and  include  both  civil  and 
military  uses.  The  growth  of  civilian  cellular  networks  makes  geolocation  ideal  for 
services  such  as  targeted  advertising  [1],  E911  [2],  roadside  assistance,  and  other 
emergency  services.  Military  applications  include  the  geolocation  of  communication 
systems,  GPS  jammers,  radars,  and  search  and  rescue  efforts. 

The  applications  which  require  geolocation  lend  themselves  well  to  the  space 
environment  for  several  reasons.  A  space-borne  geolocation  system  has  a  very  large  field 
of  view,  and  can  receive  signals  from  a  large  region  of  the  Earth.  There  are  no  restrictions 
on  over-flying  an  area  from  space,  as  there  would  be  from  the  air.  Satellites  do  not  require 
refueling,  and  can  stay  in  orbit  for  years  performing  the  mission  with  minimal  human 
intervention.  Surveillance  of  remote  areas  may  also  be  easier  for  satellites,  as  weather 
does  not  become  an  issue  in  space. 
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1.1  Problem  Statement 


The  goal  of  this  thesis  is  determine  if  emerging  methods  of  geolocation  (direct 
position  determination,  instantaneous  received  frequency,  time  difference  of  arrival,  angle 
of  arrival)  can  be  used  in  a  nanosatellite  based  system.  The  performance  these  geolocation 
methods  is  to  be  compared  by  modeling  and  simulating  the  space  environment.  A  set  of 
scenarios  is  to  be  constructed  to  provide  different  circumstances  in  which  the  modeled 
satellites  would  perform  the  geolocation  mission. 

1.2  Current  Research 

A  set  of  satellites  which  receive  a  signal  from  the  ground  will  have  several 
parameters  which  can  be  measured  to  estimate  the  location  of  the  transmitter.  The  time  at 
which  the  signal  is  received  by  the  satellite  can  be  used  in  a  time  of  arrival  (TOA) 
algorithm,  and  the  time  differential  between  signal  reception  can  be  used  to  determine  the 
time  difference  of  arrival  (TDOA).  The  phase  difference  in  an  array  of  antennas  can  be 
used  to  determine  the  angle  of  arrival  (AOA).  The  frequency  difference  of  arrival  (FDOA) 
can  be  measured  to  estimate  the  coordinates  of  the  emitter.  All  of  these  parameters  can  be 
measured  and  then  used  to  estimate  the  position  of  the  emitter.  This  is  a  two-step 
approach,  which  is  the  classical  method  for  geolocating  a  transmitter. 

Recently,  single-step  methods  have  been  suggested  to  improve  performance  of  the 
geolocating  system.  Single-step  methods  gather  the  sampled  signals  at  a  single  node  or  a 
separate  processor  to  estimate  the  position  of  the  transmitter  without  extracting  the 
parameters  separately.  The  Direct  Position  Determination  (DPD)  method  is  a  promising 
new  area  of  research  that  provides  improved  performance  in  the  presence  of  noise. 

Research  into  methods  which  are  able  to  determine  a  transmitter’s  position  with  a 
single  moving  sensor  has  gained  interest  recently  as  well.  The  Instantaneous  Received 
Frequency  (IRF)  method  is  a  new  technique  which  uses  the  Doppler  of  a  moving  sensor  to 
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estimate  the  position  of  the  transmitter.  The  method  requires  only  a  single  satellite  with  a 
single  antenna,  which  is  ideal  for  the  small  space  available  in  a  Cubesat. 

1.3  Scope  and  Application 

This  research  investigates  several  possible  methods  for  implementing  a  geolocation 
system  in  a  set  of  nanosatellites.  Several  geolocation  methods  are  discussed,  and 
simulations  of  their  use  in  the  space  environment  are  presented.  This  research  seeks  to 
provide  baseline  simulations  that  compare  the  performance  of  several  geolocation 
methods.  A  set  of  variables  will  be  included  to  assess  the  usefulness  of  the  algorithms  in 
differing  conditions.  This  research  is  performed  to  inform  a  nanosatellite  designer  on  the 
effectiveness  of  these  algorithms  if  a  geolocation  mission  is  desired. 

1.4  Assumptions 

There  are  some  limitations  that  are  present  in  this  research.  The  transmitter  is 
assumed  to  be  stationary,  and  its  location  is  generally  known,  within  a  bounding  box  that 
is  definable  in  the  simulation.  For  example,  the  location  may  be  assumed  to  be  in  or  near  a 
city,  so  the  bounding  box  could  be  modeled  as  1000  km^.  The  signal  structure  is 
unknown,  but  the  transmitted  frequency  is  known.  It  is  assumed  that  there  are  no  other 
transmitters  within  the  target  box  which  are  transmitting  on  the  same  frequency,  or 
alternately  that  the  geolocation  system  is  able  to  distinguish  the  signal  of  interest  from 
others  received.  The  signal  is  assumed  to  reach  the  satellite  with  enough  power  to  be 
received  and  processed  effectively. 

1.5  Methodology 

A  number  of  current  geolocation  methods  are  modeled  to  simulate  the  reception  and 
processing  of  signals  on-board  a  satellite.  The  simulated  signals  are  attenuated,  delayed, 
and  frequency-shifted  to  represent  a  real-world  signal  passing  through  the  atmosphere 
accurately.  Once  the  signals  are  simulated  as  being  received  by  the  sensors,  they  are 
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processed  using  one  or  more  of  several  possible  geolocation  algorithms.  These  methods 
are  all  presented  in  current  literature  on  this  subject,  though  use  of  these  algorithms  in  a 
satellite  or  at  the  standoff  range  of  a  satellite  may  be  a  new  application.  The  methods  are 
simulated  in  Matlab  to  provide  an  estimate  of  the  transmitter’s  position,  which  is 
compared  to  the  true  position  to  find  the  error  and  variance  of  the  method  over  multiple 
runs.  In  addition  to  the  average  error  of  the  methods,  other  factors  such  as  performance  at 
different  signal  to  noise  ratios  (SNR’s),  satellite  constellation  configurations,  and 
performance  factors  are  explored. 

1.6  Research  Benefits 

There  are  currently  no  nanosatellites  on  orbit  which  are  performing  a  geolocation 
mission.  This  research  intends  to  take  a  step  towards  enabling  such  a  mission  to  exist.  The 
algorithms  presented  in  subsequent  chapters  of  this  thesis  can  be  put  to  use  in  a  Cubesat, 
and  can  be  further  modified  to  suit  the  needs  of  the  user.  United  States  Air  Force  (USAF) 
and  Air  Force  Institute  of  Technology  (AFIT)  researchers  can  continue  to  improve  on  the 
work  presented  herein  to  make  it  more  robust.  The  capability  to  passively  geolocate  a 
terrestrial  transmitter  from  space  will  aid  in  mapping  out  Command,  Control,  and 
Communications  (C3)  nodes  and  creating  Operational  Plans  (OPLANs).  Geolocation  of 
hostile  jamming  sites  will  aid  the  operational  community  by  promoting  rapid 
reconfiguration  of  navigation  aids.  Search  and  rescue  operations  may  be  enhanced, 
especially  in  remote  areas. 

1.7  Thesis  Organization 

Chapter  I  introduces  the  topic  of  geolocation,  and  provides  an  overview  of  the 
research  conducted.  Chapter  II  provides  background  information  and  addresses  the 
current  state  of  geolocation  technology.  Chapter  III  introduces  the  methodology  used  to 
simulate  the  satellite  working  to  geolocate  the  signal  source  in  the  expected  environment. 
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Chapter  IV  discusses  the  results  of  simulations.  Chapter  V  presents  conclusions  from  the 
research  and  details  possible  avenues  of  future  research. 
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II.  Background 


Geolocation  of  radio  frequency  (RF)  signals  can  be  defined  as  “the  problem 

of  precise  localization  (or  geolocation)  of  spatially  separated  sources  emitting 
electromagnetic  energy  in  the  form  of  radio  signals  within  a  certain  frequency  bandwidth 
by  observing  the  received  signals  at  spatially  separated  sensors”  [3].  Geolocation  of 
signals  is  an  important  capability  for  military  and  civilian  users  alike.  Increasing  the 
accuracy  of  location  estimates  of  a  target  RF  emission  will  lead  to  more  accurate 
targeting,  better  intelligence  for  the  warfighter,  and  enhanced  search  and  rescue 
effectiveness.  Civilian  programs  like  E91 1  [4]  have  already  made  use  of  geolocation 
methods  to  make  tracking  cellphones  a  reality.  The  E91 1  program  has  put  in  place 
infrastructure  which  can  calculate  a  geolocation  estimate  every  time  a  cell  phone  connects 
to  a  91 1  operator.  The  program  has  opened  the  doors  to  research  into  improving  solutions 
to  geolocation  problems.  This  research  aims  to  add  a  robust  solution  making  use  of 
nanosatellites  (“Cubesats”)  to  the  rapidly  growing  research  area. 

This  chapter  details  the  state  of  the  art  in  geolocation  techniques.  Eiterature  in  the 
field  of  geolocation  describes  the  various  methods  to  locate  the  origin  of  a  signal  both 
passively  and  actively.  This  research  will  focus  on  the  passive  geolocation  problem,  with 
emphasis  on  the  implementation  of  these  methods  in  a  geolocation  system  involving  a 
constellation  of  Cubesats.  Sections  2. 1  and  2.2  give  an  introduction  to  the  algorithms  used 
for  geolocation.  Section  2.3  details  the  mathematical  foundation  of  geolocation  solutions. 
Section  2.4  discusses  the  primary  sources  of  error  in  a  geolocation  system.  Section  2.5 
provides  some  performance  metrics  which  are  used  to  evaluate  the  algorithms.  Sections 
2.6  to  2.7  give  the  reader  a  primer  on  Cubesat  systems  and  their  capabilities. 
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2.1  Classical  Methods  for  Locating  Emitters 

There  are  a  number  of  methods  for  estimating  the  position  of  a  transmitter  given  the 
reeeption  of  its  signal  aeross  a  set  of  sensors.  There  are  four  methods  referred  to  in  this 
thesis  as  the  ‘elassieal  methods’  -  angle  of  arrival,  time  of  arrival,  time  differenee  of 
arrival,  and  the  frequeney  differenee  of  arrival  method  [4].  These  methods  take  a  two  step 
approaeh  to  formulate  a  position  estimate.  Eaeh  method  first  ealeulates  some  parameter 
based  on  the  reeeived  signal  features,  then  eomputes  a  best  fit  estimate  of  the  transmitter’s 
position  based  on  eaeh  sensor’s  reported  parameters.  Classieal  methods  are  diseussed 
further  in  this  seetion,  and  Seetion  2.2  gives  more  information  on  alternative  methods. 
Subseetion  2.1.1  presents  the  AOA  method,  2.1.2  diseusses  the  TOA  method,  2.1.3 
introduees  the  TDOA  method,  and  2.1.4  presents  the  FDOA  method. 

2.1.1  Angle  of  Arrival. 

The  AOA  method  solves  the  geoloeation  problem  by  determining  a  Line  of  Bearing 
(LOB)  from  eaeh  sensor  to  the  target  [5].  The  LOB  ean  be  ereated  by  measuring  the  phase 
differenees  between  a  set  of  reeeived  signals  from  a  phased  array  antenna.  A  steering 
veetor  defines  the  sensor’s  best  guess  at  the  preeise  angle  at  whieh  the  signal  arrived.  The 
method  requires  at  least  two  LOBs  for  a  position  estimate  to  be  made.  Alternately,  a 
single  mobile  reeeiver  ean  be  used.  A  mobile  reeeiver  has  the  ability  to  simulate  multiple 
stationary  reeeivers  by  taking  measurements  in  different  positions  over  time,  provided  the 
position  and  veloeity  of  the  sensor  at  eaeh  time  epoeh  are  known.  The  LOBs  interseet  at 
the  estimated  point  of  signal  origin  on  a  2D  plane  or  3D  surfaee,  as  shown  in  Figure  2.1. 
Additionally,  the  sensors  for  a  satellite -based  system  are  expeeted  to  be  several  hundred 
kilometers  away  from  the  transmitter  at  their  elosest  point.  Longer  distanees  between 
sensor  and  transmitter  will  reduee  the  error  margin  that  ean  be  allowed  in  the  angle 
measurement.  The  AOA  method  requires  a  minimum  of  one  satellite  and  does  not  require 

time  synehronization  between  the  reeeiver  and  the  transmitter  for  an  aeeurate 
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B 

Figure  2.1:  The  Angle  of  Arrival  method  requires  at  least  two  lines  of  bearing  to  calculate 
the  point  of  intersection  [4] . 


measurement.  Both  of  these  characteristics  make  it  an  appealing  method  for  an 
inexpensive  passive  geolocation  system.  One  disadvantage  of  the  AOA  method  is  that 
position  estimates  are  significantly  degraded  when  the  line  of  sight  (LOS)  signals  are 
interfered  with  by  the  reception  of  Non-Line  of  Sight  (NLOS)  signals.  In  the  context  of 
this  research,  NLOS  is  not  believed  to  be  a  vital  concern  because  the  detection  will 
generally  be  done  by  space-based  receivers  with  direct  line  of  sight  to  the  transmitter.  The 
AOA  method  is  studied  primarily  for  its  ability  to  geolocate  with  a  single  satellite. 

2.1.2  Time  of  Arrival. 

Another  method  to  geolocate  an  RF  emitter  is  TOA,  which  measures  the  time  a  signal 
arrives  at  a  receiver  [4].  Given  the  signal  travels  at  the  speed  of  light,  c,  the  distance 
traveled  can  be  calculated  from  the  propagation  time.  A  circle  of  infinitely  many  possible 
emission  locations  can  be  drawn  around  the  receiving  sensor  based  on  the  calculated 
distance  that  the  signal  has  traveled.  Adding  receivers  will  reduce  the  ambiguities  of  the 
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method  by  creating  intersections  along  the  circle.  Thus,  a  group  of  three  receivers  can 
pinpoint  an  emitter  in  2-D  as  shown  in  Figure  2.2. 


Figure  2.2:  The  Time  of  Arrival  method  uses  a  signal’s  propagation  time  to  measure  the 
distance  the  signal  travelled.  The  position  estimate  is  the  point  in  space  where  the  three 
location  estimate  rings  meet  [4]. 


The  TOA  method  requires  more  receivers  than  the  AOA  method  -  three  for  a 
two-dimensional  position  estimate,  and  four  receivers  for  a  three-dimensional  estimate. 
Additionally,  the  TOA  method  requires  a  detailed  knowledge  of  the  signal  structure  to 
estimate  the  signal  propagation  time.  TOA  depends  on  precision  timing  between  the  signal 
source  and  the  receiver.  This  research  focuses  on  the  geolocation  of  non-cooperative 
receivers,  which  requires  the  assumption  that  accurate  timing  between  transmitter  and 
sensor  is  unavailable.  For  this  reason  the  TOA  method  will  not  be  used  in  this  research. 
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2.1.3  Time  Difference  of  Arrival. 

The  TDOA  method  is  perhaps  the  most  widely  used  method  for  geolocation  [6],  [7]. 
This  method  uses  the  elapsed  time  between  signal  receptions  at  two  or  more  receivers  to 
estimate  the  emitter’s  position.  Rather  than  drawing  location  estimate  rings  around  a 
receiver  (as  in  the  TOA  method),  the  TDOA  method  creates  hyperbolic  curves  to  represent 
the  possible  locations.  The  equation  for  the  hyperboloid  R,j  is  as  follows  [8]: 

Ri,j=  yj  (xi  -  Xtf  +  (y,-  -  ytf  +  {zi  -  Ztf  -  yj(xj  -  Xtf  +  (y,  -  y,f  +  {zj  -  Ztf  (2.1) 

where  the  points  {x^yuzd  and  (xj,yj,Zj)  are  the  coordinates  of  the  receivers,  i  and  j,  and 
(xt,yt,Zt)  are  the  coordinates  of  the  transmitter.  Figure  2.3  illustrates  an  example  of  how 
these  hyperboloids  can  be  used  for  geolocation.  Like  the  TOA  method,  TDOA  requires 


TDOAbq 


Figure  2.3:  The  Time  Difference  of  Arrival  method  determines  location  by  measuring  the 
difference  in  reception  time  between  two  sensors,  creating  hyperbolas  of  possible  location 
estimates  [4]. 


three  sensors  for  a  two-dimensional  estimate  and  four  sensors  for  a  three-dimensional 
estimate.  The  intersection  of  two  or  more  hyperboloids  gives  the  estimate  of  the  signal 
source. 
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TDOA  has  many  benefits  for  use  in  passive,  space-based  geolocation  systems  in 
comparison  with  the  other  methods  discussed.  TDOA  does  not  require  a  common  clock 
between  the  transmitter  and  receiver,  making  transmitter  timing  inconsequential.  The  only 
relevant  time  measurements  are  those  between  receivers,  so  it  is  critical  that  receiver 
timing  remains  synchronized  [9].  Another  benefit  of  TDOA  is  the  simplification  of  the 
LOS/NLOS  relationship.  So  long  as  there  is  one  line  of  sight  measurement,  all  others  will 
be  neglected  if  the  leading  edge  of  the  waveform  is  used  to  determine  the  TOA.  However, 
problems  arise  when  there  is  no  LOS.  Position  estimation  is  severely  degraded  when  no 
direct  path  exists  because  the  receiver  will  treat  the  first  signal  it  receives  as  the  LOS 
signal.  Without  LOS,  the  first  signal  will  originate  from  a  multipath  reflection.  With  a 
space-based  system,  NLOS  scenarios  are  assumed  to  be  minimized,  but  can  still  occur  in 
urban  canyons.  Various  NLOS  mitigation  techniques  and  algorithms  are  presented  in  [4], 
and  NLOS  mitigation  techniques  are  not  a  consideration  in  this  research.  TDOA  is  studied 
for  its  well  developed  algorithms  and  common  use  in  systems  with  multiple  sensors  at 
differing  baseline  distances. 

2.1.4  Frequency  Difference  of  Arrival. 

The  FDOA  method,  also  referred  to  as  Differential  Doppler  (DD),  measures  the 
difference  in  frequency  of  a  received  signal  between  two  or  more  receivers,  much  like  the 
TDOA  method  does  with  time  [4].  FDOA  methods  therefore  require  movement  of  either 
the  signal  source  or  the  receivers.  The  shift  in  frequency  is  used  to  obtain  a  distance 
estimate  by  using  the  relation  between  the  signal’s  speed,  c,  and  the  length  of  the 
waveform.  FDOA  estimation  can  be  visualized  by  drawing  wavelength  curves  around  the 
receivers  and  finding  points  of  intersection,  as  shown  in  Figure  2.4. 

Noise  will  degrade  the  accuracy  of  the  estimate.  The  velocity  of  the  emitter  and/or 
receivers  can  also  cause  poor  performance  if  the  movement  is  not  at  a  sufficient  speed  to 
generate  a  measurable  Doppler  shift.  Velocity  issues  will  not  be  a  significant  factor  in  the 
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Figure  2.4:  The  Frequency  Difference  of  Arrival  method  requires  a  moving  transmitter  or 
receiver  to  make  use  of  the  Doppler  shift  in  the  received  signal  [10]. 


research  in  this  thesis,  as  a  satellite  in  LEO  will  travel  at  approximately  7,500  km/s, 
ensuring  a  measurable  difference  in  estimated  frequency  over  time.  In  many  systems, 
FDOA  techniques  are  included  as  a  supplement  to  methods  such  as  TDOA  and  AOA, 
rather  than  as  a  stand-alone  solution  [7]. 

2.2  Alternative  Methods 

In  addition  to  the  classical  methods  discussed  in  the  previous  sections,  new  methods 
have  been  developed  to  improve  on  some  of  the  inefficiencies  and  pitfalls  of  the  classical 
methods.  One  such  method  is  the  Direct  Position  Determination  (DPD)  method,  which 
uses  a  single  step  cost  function  approach  to  process  the  signals  and  determine  an  estimate. 
DPD  and  other  alternative  methods  are  often  founded  on  classical  methods,  but  have 
performance  enhancements  that  distinguish  them.  Other  alternative  methods  abandon  the 
classical  approach  altogether  in  favor  of  a  new  paradigm. 

A  common  approach  is  to  use  a  hybridization  of  the  classical  methods  to  improve 
accuracy  and  robustness.  By  combining  TDOA  and  AOA  methods,  for  example,  a  system 
can  generate  estimates  based  on  timing  and  lines  of  bearing,  then  use  a  best  fit  model  to 
combine  the  two  sets  of  measurements.  Such  a  system  will  be  more  robust  by  taking 
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advantage  of  the  benefits  of  both  AOA  and  TDOA,  while  mitigating  the  negative  aspeets 
of  eaeh.  In  addition  to  hybrid  systems,  other  methods  have  been  developed  to  minimize 
the  speeifie  disadvantages  of  some  systems  in  order  to  improve  the  results  [11]  -  [14].  The 
Direet  Position  Determination  method  is  a  new  approaeh  that  is  rising  in  popularity  [15]  - 
[18]  due  to  its  single  step  solution  and  improved  performanee  at  low  SNR,  relative  to  the 
elassieal  approaehes. 

2.2.1  Direct  Position  Determination. 

The  DPD  method  improves  upon  the  elassieal  methods  through  eentralized 
proeessing  of  the  reeeived  signal.  Weiss  posits  that  there  is  error  inherent  in  elassieal 
algorithms  beeause  reeeived  signals  are  not  related  to  one  another  in  proeessing  [18]. 
Signals  are  proeessed  separately,  then  eompared  to  find  their  interseetion.  In  order  to 
eliminate  this  error,  the  signals  must  be  refereneed  to  one  sensor.  The  DPD  method 
improves  position  aeeuraey,  most  signifieantly  at  low  SNR,  but  at  the  eost  of  inereased 
eomputational  eomplexity  at  a  single  node.  Computational  eomplexity  presents  problems 
in  implementation  due  to  the  inereased  eost  for  required  proeessing  power,  as  well  as  the 
relianee  on  a  single  sensor  or  ground  node  for  the  majority  of  the  work.  Pourhomayoun 
and  Fowler  have  modified  the  proeess  to  distribute  the  ealeulations  among  all  sensors  in 
the  solution  [16].  By  using  a  Gershgorin  Dise,  eigenvalues  ean  be  approximated  to  reduee 
the  ealeulations  required.  Additionally,  they  propose  an  algorithm  that  approximates  the 
DPD  while  using  a  distributed  proeessing  model  to  aeeomplish  the  ealeulations.  Their 
results  show  that  the  method  is  nearly  equivalent  to  the  original  DPD  method,  and  still 
more  aeeurate  than  elassieal  solutions. 

The  DPD  method  has  been  well  developed  for  oases  when  the  distanoe  between  the 
transmitter  and  reoeivers  is  less  than  10  km.  Generally,  the  solution  spaoe  is  in  two 
dimensions.  The  applioation  of  the  DPD  method  to  a  spaoe  system  whioh  takes 
measurements  as  far  away  as  1000  km,  and  with  a  solution  in  three  dimensions  has  not 
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been  identified  in  literature  up  to  this  point.  This  thesis  is  believed  to  present  the  first 
performance  consideration  of  this  expanded  scenario. 

2.2.2  Instantaneous  Received  Frequency. 

The  instantaneous  received  frequency  (IRF)  method  is  related  to  the  FDOA  method 
described  in  Section  2.1.4,  but  instead  of  computing  the  difference  in  Doppler  frequencies 
of  a  signal  received  by  different  sensors,  the  instantaneous  received  frequency  is 
calculated.  A  moving  sensor  with  a  stationary  transmitter  will  have  a  velocity  which 
changes  relative  to  the  transmitter’s  position.  The  relative  velocity  of  the  satellite  will 
change  the  frequency  of  the  received  signal.  As  the  satellite  moves  toward  the  emitter,  the 
received  frequency  will  be  higher  than  the  transmitted  frequency.  This  phenomenon 
reverses  as  the  satellite  moves  away  from  the  transmitter,  and  the  received  frequency  will 
be  lower.  Geolocation  can  be  performed  using  this  information  by  matching  the  measured 
instantaneous  frequency  to  Doppler  values  calculated  from  a  grid  of  possible  emitter 
locations.  Nelson  and  McMahon  [19]  present  a  method  which  makes  use  of  these 
principles.  The  method  is  attractive  for  a  nanosatellite  mission  because  it  can  work  with  a 
single  receiver  equipped  with  a  single  payload  antenna.  This  research  presents  the  first 
simulated  performance  of  the  IRF  method  in  a  scenario  with  space-borne  sensor  platforms. 

2.3  Geolocation  Problem  Formation  and  Solution 

Whether  using  a  classical  two-step  approach  or  a  single  step  approach,  the 
parameters  needed  to  solve  for  the  transmitter’s  position  remain  embedded  in  the  signal 
received  at  the  sensors.  Torrieri  [5]  details  an  iterative  solution,  and  Ho  and  Chan  [8] 
provide  a  closed  form  solution  for  the  classical  methods.  Weiss  [18]  developed  the  single 
step  method  which  will  be  used  for  the  portion  of  this  research  dealing  with  DPD.  Nelson 
and  McMahon  [19]  provide  the  mathematical  background  for  the  IRF  method  which  is 
utilized  in  this  research  as  well.  The  following  sections  will  delve  into  the  mathematical 
basis  for  each  of  these  methods  in  more  detail. 
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2.3.1  Classical  Step  1:  Parameter  Calculation. 

The  classical  solutions  calculate  a  parameter  such  as  the  range  or  the  angle  between 
the  receiver  and  the  transmitter.  Torrieri  divides  the  estimation  methods  into  two 
categories:  hyperbolic  location  systems  and  direction-finding  location  systems  [5]. 
Hyperbolic  systems  use  methods  such  as  TDOA  and  FDOA  to  locate  a  transmitter  by 
processing  signal  arrival  time  measurements  at  three  or  more  locations.  Direction  finding 
systems  use  methods  that  create  bearing  lines  either  at  multiple  stations  or  using  the  same 
station  at  different  locations.  A  mobile  receiver  with  a  known  trajectory  is  presumed  for 
single  receiver  calculations. 

2.3.1. 1  Hyperbolic  Systems. 

Hyperbolic  location  systems  compare  the  timing  of  signals  received  at  spatially 
separated  sensors  to  create  a  hyperbola  of  possible  emission  locations  with  focii  at  the 
sensors  [5].  Adding  stations  determines  the  position  estimate  by  finding  where  hyperbolas 
cross.  An  important  aspect  of  this  method  is  that  the  need  for  finding  the  transmit  time,  to, 
is  eliminated  because  the  only  timing  that  is  needed  is  the  reception  time  difference 
between  the  two  sensors,  r.  One  exception  is  if  the  transmitter  is  producing  pulses  rather 
than  a  continuous  waveform  -  in  that  case,  care  must  be  taken  to  ensure  that  the  time 
between  measurements  is  smaller  than  the  time  between  pulses.  Equations  to  describe  the 
signals  that  two  sensors  receive,  xi(t)  and  X2(t),  are  as  follows  [5]: 


Xi(t)  =  s(t)  +  Wi(t) 

(2.2) 

X2(t)  =  as(t  -  r)  -1-  W2(t) 

(2.3) 

A  signal  s(t)  arrives  with  a  time  shift  r  at  one  sensor,  an  attenuation  factor  a,  and  some 
noise  w  which  will  be  modeled  as  white  Gaussian  noise.  The  value  of  t  can  be  measured 
by  comparing  the  reception  time  of  one  sensor  with  another.  Therefore  (2.2)  and  (2.3)  can 
be  replaced  by  establishing  the  time  difference  between  sensors  i  and  /  -l-  1  as  -  t;_i.  The 
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resulting  equation  is  the  following  [5]: 

AD 

ti  -  ti+i  = - +  w  (2.4) 

c 

where  AD  =  Di  -  D,  is  the  distance  from  receiver  i  to  the  transmitter,  and  w  is  noise. 
In  (2.4)  the  time  values  can  be  measured.  The  noise  value  can  be  modeled  stochastically, 
and  is  assumed  to  be  white  Gaussian  noise.  The  only  unknown  in  (2.4)  is  AD.  By 
measuring  the  reception  time  at  each  receiver,  solving  for  AD  is  straightforward,  and  the 
method  is  called  the  leading  edge  method  [9].  The  leading  edge  method  is  very 
computationally  inexpensive,  but  neglects  many  possible  error  scenarios.  One  would  have 
to  assume  that  both  sensors  receive  the  same  signal,  without  any  difference  in  atmospheric 
conditions,  and  with  perfect  line  of  sight.  Such  conditions  are  generally  not  the  case  in 
mission  scenarios,  so  more  complex  methods  are  employed. 

A  more  robust  method  is  to  find  the  peak  of  the  cross  correlation  between  two 
received  signals.  The  cross  correlation  measures  the  similarity  between  two  signals  as  a 
function  of  time.  This  cross-correlation  enables  estimation  of  the  time  difference  between 
signals.  The  equation  for  determining  the  time  delay  r^eiay  using  the  cross  correlation  is  as 
follows  [20]: 

Tde/ay  =  aTg  max  [(rj  ★  r2)(0]  (2.5) 

t 

where  vi  and  r2  are  the  two  received  signals,  and  the  *  operator  represents  the  cross 
correlation  operation.  This  equation  assumes  a  continuous  signal,  but  the  process  is 
generally  assumed  to  take  place  with  sampled  signals.  A  sampled  signal  is  one  which  is 
digitized  by  taking  discrete  samples  of  a  continuous  signal  such  that  the  characteristics  of 
the  original  signal  are  recoverable.  An  accurate  result  requires  the  use  of  time  steps 
between  samples  that  are  small  enough  to  achieve  the  resolution  required  for  an  accurate 
representation  of  the  original  signal.  The  required  resolution  for  localization  purposes 
depends  on  the  distance  between  the  sensors,  the  range  to  the  target,  and  the  bandwidth  of 


16 


the  signal  of  interest  (SOI).  Resolution  is  especially  important  in  space  based  systems 
because  the  range  is  so  long. 

In  practice,  to  obtain  the  time  difference  of  arrival  through  cross  correlation,  a  Fourier 
transform  is  performed  on  each  of  the  received  signals,  then  the  conjugate  of  one  signal  is 
multiplied  with  the  other  [21].  By  doing  so,  the  cross  spectral  density  (CSD)  of  the  two 
signals  is  obtained.  The  CSD  is  then  normalized  and  the  inverse  Fourier  transform  is 
computed.  The  resulting  vector  will  have  a  maximum  at  the  sample  associated  with  the 
time  difference  of  arrival  between  the  two  sensors.  The  sampling  rate  At  is  known,  so  the 
TDOA  can  be  calculated  using  the  difference  in  sample  numbers.  Normal  practice 
involves  finding  the  TDOA  associated  with  a  common  sensor  -  i.e.,  have  one  sensor  serve 
as  the  basis  for  all  TDOA  measurements  across  the  network.  The  method  employed 
provides  N  -  I  TDOA  measurements. 

Another  way  to  establish  the  cross  correlation  between  two  signals  is  the  Cross 
Ambiguity  Function  (CAF)  [21].  The  CAF  is  the  method  used  in  this  research  because  the 
DPD  algorithm  relies  upon  this  calculation.  Avoiding  the  duplication  of  effort  in  this  case 
reduces  computation  time.  The  peak  magnitude  of  the  CAF  matrix  is  a  maximum 
likelihood  estimate  for  the  TDOA  [16].  Stein  provides  the  following  equation  for  the  CAF 
[21]: 


RirJd)  = 


^  f  r\{f  -  f,)r2{f)e-R^f^df 

^  J 


(2.6) 


where  T  is  the  collection  time  span,  fd  is  the  Doppler  shift,  r  is  the  time  delay,  and  the  * 
represents  the  complex  conjugate. 

Other  approaches  are  available  and  well  documented,  and  typically  involve  a  grid 
search  or  a  Taylor  Series  expansion  [5].  These  approaches  may  require  less  processing, 
but  are  not  guaranteed  to  converge,  and  may  require  some  initialization. 
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23.1.2  Direction  Finding  Systems. 

Multiple  bearing  measurements  of  passive  direction-finding  sensors  can  be  combined 
by  a  direction  finding  system  to  produce  an  estimate  of  transmitter  position.  Once 
directions  are  determined  from  at  least  two  sensors,  the  LOBs  are  propagated  and  the 
solution  is  indicated  by  the  point  in  space  where  these  two  (or  more)  lines  cross.  As  with 
hyperbolic  location  systems,  noise  limits  the  accuracy  of  this  method.  Two  dimensional 
position  estimation  is  often  called  triangulation.  The  following  equation  defines  how  to 
determine  a  LOB  estimate,  0,  from  sensor  to  transmitter,  neglecting  noise  [5]: 


6i  =  cos 


-1 


Xt  -  Xi 


,  0  <  9i  <  n 


(2.7) 


[  ^(xt  -  Xif  -r  (jt  -  ytf  {zt  -  ZiY  \ 
where  \_Xt,yt,ZtY  are  the  transmitter’s  coordinates  and  are  the  sensor’s 

coordinates.  It’s  apparent  that  0,  has  an  ambiguity  such  that  the  quadrant  cannot  be  derived 
outright.  Knowledge  of  the  environment  and  the  target’s  approximate  location  can  be 
applied  to  eliminate  the  ambiguity.  For  example,  a  satellite  looking  at  Earth  can  eliminate 
two  quadrants  altogether  by  focusing  on  the  two  quadrants  where  the  Earth  resides. 

2.3. 1.3  Phase  Interferometry. 

Phase  interferometry  is  a  method  for  determining  the  angle  of  arrival  by  measuring 
the  phase  difference  of  an  incoming  signal  between  two  or  more  antennas  [22].  If  both  the 
distance  between  antennas  and  the  carrier  frequency  of  the  incoming  signal  are  known,  the 
phase  relationship  of  the  received  signals  can  be  determined  and  the  angle  of  arrival 
becomes  apparent.  Eigure  2.5  shows  the  general  principle  of  phase  interferometry.  This 
method  may  provide  very  accurate  results  with  well  tuned  hardware  and  appropriate 
signal  sampling  frequencies.  However,  the  accuracy  is  not  without  some  drawbacks. 

Phase  interferometry  requires  at  least  two  antennas,  and  this  will  only  provide  about  120 
degrees  field  of  view  [22].  In  order  to  attain  a  full  360  degree  field  of  view,  along  with 
both  azimuth  and  elevation  information,  a  minimum  of  four  antennas  are  required  [22]. 

An  additional  problem  with  using  phase  interferometry  exclusively  results  from  receiving 
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multiple  signals  arriving  from  different  angles.  In  order  to  mitigate  some  of  the  issues  with 
using  phase  interferometry,  the  MUSIC  algorithm  is  employed  for  this  researeh. 


Figure  2.5:  Phase  Interferometry  measures  the  phase  differenee  between  multiple  antennas 
when  a  plane  wave  hits  the  array. 


2.3. 1.4  The  MUSIC  Algorithm. 

The  MUSIC  algorithm  is  a  oornmon  teehnique  in  direetion  finding  applieations  [23]  - 
[25].  The  teehnique  uses  the  eovarianee  matrix  obtained  from  M  signals  reeeived  by  an 
antenna  array  to  determine  the  angle  at  whieh  the  signals  arrived.  The  basie  signal  model, 
x(t),  is  as  follows  [26]: 

x(t)  =  As(0  +  w(t)  (2.8) 

where  x(t)  is  of  length  L,  equal  to  the  number  of  reeeiving  antennas.  A  is  a  veetor  of  M 
steering  veetors,  a(0),  s(t)  is  an  M  x  1  veetor  of  the  souree  signals,  and  w(t)  is  an  L  x  1 
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noise  vector.  The  steering  vector  a(9)  is  defined  as  follows  [26] : 


1 


m  = 


-ja>- 


(2.9) 


.  (L-l)rfsin6^ 

e  ^ 


where  cjj  is  the  frequency,  d  is  the  distance  between  antenna  array  elements,  and  c  is  the 
speed  of  light.  In  order  to  use  (2.8),  the  correlation  matrix  must  be  defined  [26]: 


Rx  =  E{x{t)x^{t)]  =  ARsA"  +  £{w(0w^(0}  (2.10) 


where  Rs  =  E  |s(0s^(0}  =  diag  |  cTj, . . . ,  cr^j,  the  variance  of  the  signal  at  each  receiver. 
The  special  case  where  the  elements  of  the  noise  vector  w  are  mean  zero,  variance  ctq, 
changes  the  value  of  E  |w(0w^(0|  to  cTqI  [26].  So  long  as  L>  M,  the  matrix  ARjA^  is 
singular.  Since  det|^ARsA^j  is  0,  cTq  must  be  an  eigenvalue  of  the  correlation  matrix  Rx.  In 
order  to  execute  this  procedure,  first  find  the  eigenvalues  and  eigenvectors  of  the  Lx  L 
correlation  matrix,  then  partition  them  into  two  sections.  The  first  are  the  eigenvectors 
associated  with  non- zero  valued  eigenvalues,  called  U^.  The  second  section  consists  of  the 
eigenvectors  associated  with  zero  eigenvalues,  called  U„.  The  eigenvectors  and  U„ 
represent  the  signal  subspace  and  the  noise  subspace,  respectively,  and  are  represented  as 
follows  [26]: 


U.  U„ 


Ml  ...  Um  Um+\  ■  ■  ■  Ul 


(2.11) 


where  M  is  the  number  of  unique  signals  that  are  being  received. 

The  steering  vector  defined  in  (2.9)  is  in  the  signal  subspace,  which  by  definition  is 
orthogonal  to  the  noise  subspace.  It  is  therefore  known  that  a^(0)Un  =  0.  If  a  search 
through  all  possible  angles  from  0°  to  90°  is  performed,  the  pseudo  spatial  spectrum  in 
which  the  correct  AOA  lies  can  be  mapped.  The  following  equation  for  the  pseudo 
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spectrum  P{6)  is  then  used  [26]: 


p{e)  = 


1 

l|a«(0)U„||2 


(2.12) 


where  U„  is  a  matrix  of  noise  subspace  eigenvectors.  Wherever  there  exists  a  peak  in  P{6), 
there  also  lies  a  potential  AOA  for  the  arriving  signal.  For  3D  estimates,  this  procedure 
must  be  completed  for  both  azimuth  and  elevation. 

2.3.2  Classical  Step  2:  Lateration/Angulation  Techniques. 

Once  the  results  from  each  sensor’s  calculations  are  complete,  they  can  be  sent  to  a 
central  node  for  further  processing.  For  most  methods,  the  goal  of  the  second  step  is  to 
find  the  places  where  lines  or  hyperbolas  cross,  representing  potential  transmitter 
locations.  Ideally,  all  lines  will  cross  at  a  single  point,  but  with  the  presence  of  noise  in 
real-world  systems,  this  is  not  the  case.  At  the  central  node,  an  estimate  for  the  emitter’s 
location  in  the  presence  of  noise  is  formed  by  way  of  a  best  fit  approach,  generally  either 
Maximum  Likelihood  Estimation  (MLE)  or  Eeast  Squares  Estimation  (ESE)  techniques. 

2.3.2. 1  Maximum  Likelihood  Estimation. 

Position  estimation  in  the  presence  of  noise  requires  the  use  of  a  statistical  model 
rather  than  the  deterministic  equation  it  would  be  in  the  absence  of  noise.  A  set  of  sensors 
which  obtains  Ns  samples  of  the  signal  will  yield  Ns  position  estimates,  each  of  which 
have  some  uncertainty  associated  with  them.  The  MEE  method  uses  the  inherently  limited 
number  of  samples  to  determine  the  parameters  which  maximize  the  parameters  of  the 
likelihood  function,  and  minimize  the  uncertainty.  Casella  and  Berger  give  the  following 
definition  of  MEE  [27]: 


Eor  each  sample  point  x,  let  6(x)  be  a  parameter  value  at  which  L(6\x)  attains 
its  maximum  value  as  a  function  of  6,  with  x  held  fixed.  A  maximum 
likelihood  estimator  (MEE)  of  the  parameter  6  based  on  a  sample  X  is  6(X). 
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The  solution  to  an  MLE  problem  is  typically  found  by  finding  the  instances  when  the 
derivative  of  the  likelihood  function  is  zero.  The  derivative  will  provide  a  maximum, 
minimum,  or  point  of  inflection.  A  second  derivative  can  provide  further  information  to 
define  the  maximum  or  minimum  of  the  function.  In  cases  when  multiple  maxima  are 
available,  such  as  when  there  are  multiple  satellites  with  differing  position  estimates,  the 
global  maximum  is  used. 

2.3.2. 2  Least  Squares  Estimation. 

The  LSE  technique  is  used  in  geolocation  systems  with  overdetermined  solutions,  i.e. 
a  solution  with  more  equations  than  unknowns.  Eor  example,  a  system  using  a 
constellation  of  five  satellites  to  estimate  a  three-dimensional  location  is  overdetermined 
because  only  four  satellites  are  required.  The  ESE  method  can  be  used  to  reduce  error  in 
the  estimate  by  computing  a  solution  that  minimizes  the  sum  of  the  squared  residual  terms 
in  a  set  of  range  estimates.  The  ESE  technique  is  commonly  used  in  AOA  methods  where 
EOBs  are  aggregated  over  a  period  of  time  or  for  an  overdetermined  set  of  sensors.  To 
complete  the  least  squares  calculations,  the  quantities  m,.  A,  and  b  are  defined  [28]: 


Vj  -  Wj 
\\Vi  -  Wu 


(2.13a) 


N 

A  =  ^  I  -  UiuJ  (2.13b) 

i=l 

N 

b  =  ^{l-UiuJ)wi  (2.13c) 

i=l 

where  v,  is  the  EOB  start  point’s  coordinates,  w,  is  the  EOB  end  point’s  coordinates, 
and  I  is  the  identity  matrix.  Using  the  definitions  in  (2.13)  the  minimum  error  solution  can 
be  determined.  The  least  squares  estimate  of  the  target’s  position  through  AOA  estimation 
is  the  solution  to  the  equation  x  =  A~^b. 

2.3.3  Closed  Form  Solution. 

Ho  and  Chan  have  developed  a  closed-form  solution  to  the  TDOA  based  geolocation 
problem,  and  address  the  issues  that  arise  with  linearization  and  iterative  methodologies 
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(such  as  those  proposed  by  Torrieri  [5]).  The  study  concerns  the  use  of  geostationary 
satellites  for  positioning  because  geostationary  satellites  are  motionless  with  respect  to  a 
point  on  the  Earth,  which  reduces  error.  The  solution  described  uses  the  geocentric 
coordinate  system.  TDOA  for  each  satellite  pair  is  defined  as  the  time  from  origin  at 
which  the  correlation  function  attains  the  largest  value.  The  solution  for  the  range  between 
the  transmitter  an  the  receiver,  r,-,  proposed  by  Torrieri  requires  linearization  using 
Taylor  Series  expansion,  shown  by  the  following  equation  [5] : 

n  =  ^J(xi  -  XtY  +  (y,-  -  ytf  +  {Zi  -  Ztf,  i  =  1,2,  ...,A  (2.14) 

where  N  is  the  number  of  satellites.  The  Taylor  Series  expansion  method  requires  a 
reasonable  initial  estimate  to  converge  (convergence  is  not  guaranteed),  and  is  typically 
processed  iteratively.  It  is  also  limited  to  critically  determined  systems,  so  extra 
measurements  are  discarded.  Ho  and  Chan’s  method  is  typically  used  with  four-sensor 
arrangements,  but  may  be  expanded  to  systems  with  only  three  sensors. 

2. 3. 3.1  Four  Sensor  TDOA  Localization. 

A  geolocation  system  using  four  sensors  is  capable  of  producing  a  three  dimensional 
estimate  of  the  target’s  position  for  each  time  instant  that  the  sensors  collectively  capture 
the  SOI.  Four  sensors  can  yield  position  computations  making  use  of  three  equations  with 
three  unknowns,  so  these  equations  are  solved  to  obtain  the  emitter’s  3D  position.  Ho  and 
Chan  [8]  have  developed  a  closed  form  method  for  finding  the  target’s  position.  Kulumani 
subsequently  streamlined  this  method  [29].  The  algorithm  is  summarized  in  the  following 
discussion. 

The  relationship  Ar,  j  between  the  time  difference  and  the  Euclidean  range  difference 
between  a  pair  of  receivers  is  given  as  follows  [8]: 

Ar;,i  =  CT,-,1  =  ri  -ri,  /  =  2,3,4  (2.15) 
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where  t,  j  is  the  time  delay  between  sensor  i  and  sensor  1,  and  r,  is  the  Euclidean  range 
from  the  target  to  the  sensor.  By  inserting  the  Euclidean  range  equation  into  (2.15),  the 
three  equations  of  interest  are  written  with  respect  to  the  common  sensor’s  position, 
UuyuZiY  [8]: 

=  ^j{x,  -  XiY  +  {yt  -  yiY  +  {zt  -  ZiY  -  V(-^f  -  +  Cv?  -  y\Y  +  (zt  -  ziY  (2.16) 

Squaring  both  sides  of  (2.16)  and  expanding  yields  the  following  equation  [29]: 

=  Y  -  2Ar,- iri  -  i  =  2, 3,4  (2.17) 

Eurther  simplification  as  well  as  putting  the  equations  into  matrix  form  and  with  respect  to 
the  variable  ri  yields  the  following  equations  [29] : 

Arj^  Ar2,i  X2x  yi.i  Z2,i  Xt  K2  Ki 

Ar^j  +2ri  Arsj  = -2  jc3,i  yaj  Z34  yt  +  Kj,  -  Ki  (2.18) 

Ar^J  [Ar4jJ  [jC4,i  y4,i  Z4,iJ  [zfj  [^4]  [^1 

where  Ki  =  /  =  1, 2, 3, 4.  The  only  unknowns  in  (2.18)  are  the  range  from  the 

transmitter  to  the  common  sensor,  ri,  and  the  transmitter’s  coordinates,  {xt  yt  ZtY ■  In  order 
to  find  the  range  ri,  the  following  equations  are  used  to  define  the  temporary  vectors  a,  (5 
and  the  values  A,  B,  and  C  [29]: 

«1  X2X  y2,l  Z2,l  Ar2,i 

«2  =  “  ys.i  Z3,1  Ar3,i  (2.19) 

«3  -^4,1  y4,l  Z4,l  Ar4i 

/?!  X2,l  y2,l  Z2,1  Ar^  J  -K2  +  K1 

Pi  ~  ~  -^3,1  y3,l  Z3,\  2  '^^3,1  “  ^3  +  ^1  (2.20) 

AJ  [-^4,1  Ki  Z4,i\  [Ar^  J  -  i^4  + 

A  =  cr j  +  q;2  +  <^3  “  1  (2.21) 
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B  =  2(ai/3i  +  02/32  +  a^/S^  -  XiOi  -  JiQ'2  -  ZiQ's) 


(2.22) 


C  =  K,-  lx, 13,  -  2y,/32  -  2z,/32  +p]+/3l+ 13]  (2.23) 

These  quantities  can  be  put  into  a  polynomial  equation  as  follows  [29] : 

0  =  Ar^  +  Br,  +  C  (2.24) 

The  roots  of  (2.24)  represent  the  possible  values  of  r,.  There  will  be  multiple  roots,  and 
the  correct  one  must  be  chosen  for  use  in  (2.18)  to  obtain  the  correct  transmitter  position. 
A  common  constraint,  which  is  used  in  this  development,  is  that  the  final  solution  must  be 
a  position  on  or  near  the  Earth’s  surface.  After  the  roots  of  (2.24)  are  input  for  the  r,  term, 
(2.18)  is  rearranged  to  solve  for  the  transmitter  position  as  follows: 

Xt  X2,,  y2,i  Z2,1  Ar2,i  J  -  K2  +  K, 

yt  =  ~  ys.i  Z3,1  Ar3,i  r,  +  \  +  K,  (2.25) 

Zt\  [x4,,  y4,i  Z4,iJ  [[Ar4jJ  [rl^-K4  +  K,\ 

The  four  satellite  TDOA  solution  proposed  by  Ho  and  Chan  [8]  and  streamlined  by 
Kulumani  [29]  can  be  shown  to  be  a  good  estimator  for  geolocation  purposes.  This 
algorithm  is  employed  for  the  research  in  this  thesis. 

2.33.2  Three  Sensor  TDOA  Localization. 

When  only  three  sensors  are  available,  the  TDOA  method  is  still  possible,  but  a 
constraint  is  required  to  determine  the  transmitter’s  position.  The  constraint  that  the 
transmitter  is  on  the  Earth’s  surface  is  again  employed,  but  with  only  three  sensors  this 
constraint  will  take  a  much  larger  role  than  in  Section  2.3.3. 1. 

The  equation  for  the  chief  constraint  is  as  follows  [29] : 

r2  =  .r?+y2  +  ^2  (2.26) 

where  {x,,  y,,  Zt\  are  the  coordinates  of  the  transmitter  in  the  Earth-centered  Earth-fixed 
(ECEE)  coordinate  frame.  Reusing  (2.17)  and  the  previously  defined  value  of  K,, 
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equations  are  placed  in  matrix  form  in  terms  of  as  follows  [29] : 


Xt  xi  yi  zi  K^+rl-  r] 

yt  =“2  -^2  Z2  K2  +  r]  -  rl  -  llyr2,iri  -  (121) 

Zt\  \xj,  yj,  Z3J  [i^3  +  r2-r^-2Ar34ri -Ar^  j 

where  [2c,y,z,]  represents  the  coordinates  of  the  satellite.  Notice  in  (2.27)  that  plays  a 
significant  role  in  the  emitter  position  solution,  in  contrast  with  (2.18),  which  does  not 
involve  at  all.  The  following  definitions  for  the  temporary  variables  a,  yS,  y,  a,y,  and  A 
through  I  aid  in  solving  for  the  range  from  the  satellite  to  the  transmitter  [29]: 

a  Ki+r] 

f3  =  K2  +  rl-Arl^  (2.28) 

yj  [K^  +  r^.-Arj^ 


an 

<3i2 

<3i3 

3Ci 

Jl 

Zl 

1 

<321 

<322 

a23  "2 

2C2 

3^2 

Z2 

<331 

<332 

<333 

-^3 

3^3 

Z3 

A  an  ai2  ai3  -1 

D  —  a2i  CI22  (^23  ~1 

G  [<231  a32  a33j  1^-1 

B  ai2  ai3  - 

-2Ar24 

E  =  a22  «23 

[-2Ar3j 

H  a32  a33 

C  Q-ii  0-12  ^13  (X 

F  =  a2l  a22  <323  /3 

I  a3i  a32  a33  y 


(2.29) 


(2.30) 


(2.31) 


(2.32) 
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The  values  defined  in  (2.28)  -  (2.32)  can  be  combined  to  form  a  quartic  polynomial  with 
respect  to  ri ,  as  follows  [29] : 


0  =  Arf  +  Br\  +  Cr^  +  Drj  +  £  (2.33) 

where 

A  =  (2.34a) 

B  =  2AB  +  2DE  +  2GH  (2.34b) 
C  =  -2;ciA  -  2yiD  -  2ziG  +  2AC  +  2DF  +  2GI  +  -  1  (2.34c) 

b  =  -2;ciB  -  2yiE  -  2ziH  +  2BC  +  2EF  +  2HI  (2.34d) 

E  =  -IxiC  -  2yiF  -  Izil  +  +  F^  +  (2.34e) 

The  roots  of  (2.33)  are  then  obtained.  There  will  be  four  roots,  which  are  the  possible 
ranges  from  the  transmitter  to  the  common  sensor.  Only  two  of  the  four  roots  will  be 
positive,  and  the  negative  roots  can  be  ignored.  The  two  positive  roots  are  substituted  back 
into  (2.27)  to  determine  the  two  possible  emitter  positions.  Additional  logic  can  be 
applied  to  determine  the  root  which  corresponds  to  the  emitter’s  position,  such  as 
comparison  of  the  solution  with  the  received  TOA’s  to  determine  the  correct  solution. 
However,  care  must  be  taken  to  ensure  that  the  root  chosen  is  on  or  near  the  surface  of  the 
Earth  at  the  given  latitude. 

2.3.4  Single  Step  Approach. 

The  DPD  algorithm  takes  a  different  approach  to  solving  for  position  with  a  set  of 
receivers.  Weiss  and  Amar  [18],  along  with  Pourhomayoun  and  Fowler  [16],  show  that 
although  each  TDOA/FDOA  estimate  is  optimal  within  each  stage  of  estimation,  the 
two-stage  approach  itself  is  not  optimal.  Error  is  introduced  by  not  relating  each 
measurement  to  a  single  sensor.  By  processing  all  of  the  signals  together  at  a  common 
node  this  error  is  eliminated. 
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The  DPD  method  transforms  the  sampled  received  signals  into  the  frequency  domain 
by  way  of  a  Fast  Fourier  Transform  (FFT)  to  separate  the  time  shift  from  the  frequency 
shift.  The  resulting  signal,  s„  is  observed  as  follows  [18]: 

s„  =  m  +  Wnitk-i)  for  6  (0, 1, ,  N,)  (2.35) 


where  s„  is  Ns  samples  of  the  received  signals  at  the  sensor,  cr„  is  a  complex  attenuation 
factor,  fd  is  the  Doppler  shift,  s  is  the  transmitted  signal,  and  w/  is  a  white  Gaussian  noise 
associated  with  the  sensor.  A  grid  of  possible  transmitter  positions  is  created  by 
mapping  the  sensor  response  to  2-D  grid  locations  in  the  matrix  Vpj  using  the  following 
equation  [16]: 


Vp,  ^  [Df  Wf sf ,  Df  Wf sf , . . . ,  D^W^s"] 


(2.36) 


where  is  the  shift  operator,  is  a  matrix  with  w(t)  values  on  the  diagonal,  and  sf  is 
a  column  vector  of  received  signal  samples.  In  order  to  find  the  position  estimate,  an 
N  xN  matrix  Qp.  =  must  be  calculated  for  each  position  on  the  grid.  An 

exhaustive  search  is  carried  out  through  the  grid  of  potential  emitter  locations.  The  point 
on  the  grid  where  Qp,  attains  a  global  maximum  eigenvalue  is  the  position  estimate. 

Pourhomayoun  and  Fowler  [16]  expand  on  this  method  by  simplifying  and 
distributing  it.  The  original  DPD  is  computationally  expensive  and  requires  one  sensor  to 
perform  the  estimation  calculations.  By  approximating  the  eigenvalue  calculations,  using 
Gershgorin’s  Theorem  [30],  the  DPD  method  becomes  much  more  attractive  for 
real-world  applications.  Stein  shows  that  a  signal  can  be  decomposed  into  its  TDOA  and 
FDOA  components  by  way  of  the  Complex  Ambiguity  Function  [31].  Pourhomayoun 
notes  that  the  Qp.  matrix  derived  in  Weiss’  work  is  merely  a  grid  of  individual  CAF 
measurements  between  sensors.  In  order  to  distribute  the  workload  of  the  algorithm,  each 
sensor  can  sample  the  received  signal  and  then  transmit  those  samples  to  the  entire  sensor 
network.  Each  sensor  then  computes  a  CAF  matrix  for  every  other  sensor  in  the  network 
and  map  the  matrices  to  the  grid  described  previously.  In  order  to  reduce  the 
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computational  load  of  computing  an  eigenvalue  for  each  point  on  the  grid,  a  number  of 
approximations  and  simplifications  can  be  made. 

First,  a  Gershgorin  Disc  is  used  to  approximate  the  eigenvalue  [16].  A  Gershgorin 
Disc  is  an  eigenvalue  approximation  which  approximates  the  eigenvalues  of  a  matrix  as 
the  diagonal  elements  of  a  row  plus  or  minus  an  uncertainty  term  defined  by  the  absolute 
value  of  the  other  elements  in  that  row  or  column,  whichever  has  a  smaller  absolute  sum. 
An  example  using  the  arbitrary  matrix  A  follows: 
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A  has  Gershgorin  Discs  D(5, 1.4),  D(3, 1.9),  D(10, 0.6),  D(-12, 1.4).  The  first  and  third 
discs  use  the  column  values,  while  the  second  and  fourth  discs  use  the  row  values. 
Gershgorin’s  Theorem  states  that  every  eigenvalue  in  matrix  A  is  located  within  one  of  the 
Gershgorin  Discs  [30].  If  the  uncertainty  term  is  sufficiently  small,  this  theorem  can  be 
used  in  determining  eigenvalues  without  the  computational  load  normally  required.  The 
eigenvalues  in  this  example  are  di  =  4.3435,  A2  =  3.5749,  /I3  =  10.0207,  A4.  =  -11.939. 

Using  Gershgorin’s  Theorem  in  conjunction  with  knowledge  of  the  grid  structure  of 
the  CAF  allows  further  simplification  [16].  The  Qp.  matrix  is  Hermitian  and  positive 
definite  [16],  so  the  columns  and  rows  are  conjugates  of  each  other.  There  is  therefore  no 
need  to  determine  whether  the  column  or  the  row  will  offer  the  smallest  uncertainty.  Now, 
all  each  sensor  must  do  is  compute  the  CAFs  in  relation  to  the  other  sensors,  find  the 
largest  value,  and  compare  it  to  the  other  sensors  for  a  final  position  estimate. 
Pourhomayoun  and  Fowler  performed  a  multitude  of  simulations  [16]  with  results  that  are 
near  the  same  quality  as  the  original  DPD  approach,  and  exhibit  less  position  estimate 
error  than  classical  methods  (TDOA)  in  low  SNR  environments. 
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2.3.5  Instantaneous  Received  Frequency  Approach. 

Nelson  and  McMahon  present  a  method  for  geolocating  an  emitter  when  its  signal  is 
observed  by  a  single  moving  receiver  [19].  The  method  begins  by  defining  the  equation 
for  the  instantaneous  observed  frequency,  [19]: 


af^\t)  =  oj  1- 


Vrif) 


(2.37) 


where  a»  is  the  transmitted  frequency  and  Vr{t)  is  the  velocity  component  between  the 
receiver  and  transmitter.  Then,  a  grid  of  potential  emitter  locations  x  is  constructed.  The 
relative  velocity  Vr  for  the  receiver  is  calculated  for  each  point  on  the  grid.  The  estimated 
instantaneous  received  frequency  x)  for  each  of  these  grid  points  is  calculated  as 
follows  [19]: 
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(2.38) 


The  authors  suggest  two  possible  ways  to  estimate  the  emitter  position  using  The 
first  method  calculates  the  estimated  IRF  for  each  grid  point,  then  determines  the  grid 
point  with  the  smallest  variance  V,  according  to  the  following  equation  [19]: 


y(jc)  =  var,  {pf^'  {t,  x)) 


(2.39) 


The  second  method  seeks  to  minimize  a  cost  function,  G  which  compares  with 
the  expected  received  frequency,  from  each  grid  point.  The  equation  for  G  is  as 
follows  [19]: 
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dt 


(2.40) 


where  G  is  the  center  frequency  of  the  signal.  Conducting  an  exhaustive  search  to 
determine  the  minimum  in  either  (2.39)  or  (2.40)  among  all  the  grid  points  provides  the 
estimated  transmitter  position  for  the  IRF  method.  Employing  an  exhaustive  search  will 
ensure  a  global  minimum  is  found. 
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2.4  Sources  of  Error  and  Mitigation 

The  accuracy  of  a  geolocation  solution  is  dependent  on  three  major  factors:  the 
accuracy  of  the  sensor’s  position,  the  accuracy  of  range/angle  estimates,  and  the  geometric 
configuration  of  the  sensors  and  the  unknown  transmitter  position  [4] .  The  first 
consideration,  the  accuracy  with  which  the  sensor’s  position  is  known,  is  directly  related  to 
the  positioning  knowledge  of  the  Cubesat.  The  precision  of  the  range  and  angle  estimates 
is  a  function  of  the  accuracy  of  the  geolocation  model  employed,  and  is  the  chief  concern 
of  the  system  design.  The  geometry  of  the  sensors  and  the  transmitter  will  be  a  concern 
due  to  position  dilution  of  precision  (PDOP).  PDOP  is  a  unitless  measure  of  a  positioning 
system’s  geometry  with  respect  to  the  range  between  the  transmitters  and  receivers.  A 
good  PDOP  is  a  low  value,  with  a  minimum  of  one,  and  a  poor  PDOP  has  a  high  value. 

2.4.1  The  Importance  of  Sensor  Position  Knowledge. 

Accurate  position  knowledge  of  a  geolocation  sensor  is  crucial  in  lowering  the  error 
of  the  geolocation  estimate.  A  distinction  should  be  made  between  position  knowledge 
and  position  control.  A  satellite’s  ability  to  know  precisely  its  position  and  attitude  is 
more  important  to  the  geolocation  solution  than  its  ability  to  control  precisely  its  position 
and  attitude.  From  500  km  away,  a  line  of  bearing  that  is  off  by  2°  will  cause  a  position 
estimate  error  of  17.5  km  [32].  This  error  can  be  compounded  by  poor  geometry  and 
timing  issues  as  well.  It  is  therefore  prudent  to  introduce  a  method  to  take  position 
knowledge  error  into  account  when  calculating  geolocation  estimates.  Ho  et  al.  provide  an 
algorithm  to  address  the  problem  of  RF  source  localization  in  the  presence  of  random 
receiver  location  errors  [33].  Ho’s  proposed  mitigation  method  has  two  stages.  First, 
linearize  the  geolocation  equations  (TDOA  is  used,  though  expansion  to  other  algorithms 
is  straightforward)  with  the  introduction  of  a  bias  parameter  related  to  the  receiver 
location.  Next,  use  the  induced  bias  to  refine  the  emitter  location  estimate  through  the  use 
of  a  weighted  least  squares  minimization  using  the  range  values  computed  in  the  first 
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stage.  Some  assumptions  are  made  along  the  way  in  the  derivation.  First,  it  is  assumed 
that  the  measurement  noise  in  the  sensor  is  independent  of  the  receiver  location  noise. 
Also,  the  algorithm  requires  knowledge  of  the  covariance  matrices  from  the  measurements 
and  from  the  receiver  locations.  These  values  would  require  some  calibration  via  in-orbit 
testing  with  a  signal  at  a  known  location.  The  result  of  the  study  was  the  derivation  of  the 
increase  in  the  Cramer-Rao  Lower  Bound  (CRLB),  and  an  algebraic  solution  to  improve 
the  estimation  accuracy  in  the  presence  of  receiver  position  knowledge  uncertainty  was 
found  to  meet  this  new  CRLB. 

2.4.2  Range  and  Angle  Estimate  Accuracy. 

The  precision  of  range  and  angle  estimation  is  largely  a  function  of  the  quality  of 
data  collected  by  the  sensors  and  the  accuracy  with  which  a  system  can  analyze  that  data. 
The  sensors  in  this  research  will  be  in  space,  450  km  above  the  Earth’s  surface.  In  order 
for  a  signal  to  reach  the  satellite’s  antenna,  it  must  have  enough  power  to  travel  that 
distance  and  still  stand  out  from  the  noise.  An  antenna  with  a  suitable  gain  pattern  is 
essential  to  detecting  the  target  signals  when  they  are  not  designed  for  communication 
with  the  sensors. 

Most  geolocation  algorithms  rely  on  an  accurate  timing  system  to  determine  the  time 
a  signal  arrives  at  a  sensor  in  order  to  calculate  either  a  distance  measurement  or  a  LOB. 
Cubesats  generally  don’t  have  the  space  available  for  a  clock  that  is  more  accurate  than  a 
crystal  oscillator,  so  the  timing  for  geolocation  related  measurements  must  be  taken  from 
larger  satellites  nearby  (most  likely  GPS  satellites).  Allan  notes  that  “...every  clock 
disagrees  with  every  other  clock  essentially  always,  and  no  clock  keeps  ideal,  or  ‘true’ 
time  in  an  abstract  sense  except  as  we  may  choose  to  define  it’’  [34].  Allan  discusses  the 
variables  which  are  prone  to  error  in  all  timing  systems,  such  as  frequency  drift,  frequency 
offset,  time  offset,  and  random  deviations.  All  of  these  variables  can  be  aggregated  into  a 
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clock  error  term,  which  is  unavoidable  in  localization  systems.  A  typical  clock  error  in 
Cubesats  that  rely  on  GPS  is  on  the  order  of  40  ns  [35]. 

2.4.3  Additional  Error  Sources. 

Additional  sources  of  error  include  multipath,  NLOS,  and  atmospheric  factors. 
Multipath  refers  to  the  propensity  of  electromagnetic  waves  to  reflect  off  of  objects 
present  in  the  path  between  the  transmitter  and  the  sensor.  Multipath  causes  multiple 
replicas  of  the  signal  to  be  received  at  the  sensor  over  time,  often  with  delayed  signals 
attenuated  due  to  the  absorption  of  energy  during  the  reflection.  The  effects  of  multipath 
include  uncertainty  in  determining  which  received  signal  is  the  original  LOS  signal,  as 
well  as  the  deterioration  of  signal  strength  and  structure.  Another  scenario  in  which 
multipath  is  an  issue  is  when  there  is  no  LOS  signal  received.  In  this  case,  there  is  an 
object  blocking  the  line  of  sight,  so  only  the  multipath  signals  arrive  at  the  sensor,  causing 
a  bias  in  the  range  estimate.  An  important  measure  of  the  effect  of  multipath  on  a  system’s 
ability  to  resolve  multipath  channels  is  the  bandwidth  [4].  As  bandwidth  is  increased,  the 
time  domain  resolution  is  increased,  which  improves  range  estimates  [4].  Narrowband 
systems  are  therefore  degraded  more  severely  than  wideband  systems.  A  model  of 
multipath  signal  h  delayed  by  r  is  as  follows  [4] : 

Lp 

Kt)  =  ^  ake^‘>”^6{t  -  Tk)  (2.41) 
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where  Lp  is  the  number  of  multipath  channels,  and  the  quantities  and  are  the 

amplitude,  phase,  and  propagation  delay  of  the  path,  respectively. 

Another  error  to  be  concerned  with  is  atmospheric  attenuation.  The  angle  at  which 
the  signal  passes  through  the  upper  atmosphere  may  have  an  effect  on  the  speed  at  which 
the  signal  travels  and  may  diffract  the  signal  as  it  travels  through  the  channel.  There  are 
many  models  that  deal  with  elements  of  atmospheric  attenuation.  Barts  [36]  and 
Panagopoulos  [37]  model  the  propagation  of  high  frequency  satellite  communication 
waves.  Wang  and  Sarabandi  [38]  examine  wave  propagation  through  foliage.  Rain  can 
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degrade  signals  in  the  K  band  as  well,  leading  to  severe  attenuation  of  the  signals  by  the 
time  they  arrive  at  the  sensors.  The  effects  of  atmospheric  attenuation  will  not  be  a  focus 
of  this  research.  Another  factor  which  affects  the  attenuation  of  a  signal  when  propagated 
through  the  atmosphere  is  the  signal’s  frequency.  A  higher  frequency  signal,  such  as  one 
in  the  X  band,  will  attenuate  differently  through  the  atmosphere  than  a  lower  frequency 
signal  such  as  one  in  the  VHF/UHF  band.  A  system’s  target  frequency  can  have  an  effect 
on  the  design  of  several  aspects  of  the  system  as  well. 

2.4.4  The  Effect  of  Frequency. 

Transmitter  frequency  is  an  important  consideration  in  a  geolocation  system.  A 
signal’s  carrier  frequency  can  limit  the  effective  range  of  communication,  the  receiver 
antenna  array  design,  the  range  between  receivers,  and  many  other  factors  which  go  into 
the  design  of  the  geolocation  system.  Airborne  receivers  typically  have  the  space  to  carry 
multiple  antenna  arrays  to  perform  a  wide  range  of  missions  with  different  frequency 
requirements.  A  small  space-borne  system  is  much  more  limited  in  both  space  available 
for  antenna  arrays  and  the  ability  to  change  hardware  with  each  new  mission.  The 
frequency  and  beamwidth  of  the  signal  may  influence  the  spacing  of  the  satellites  as  well, 
and  would  then  affect  the  PDOP.  Sensors  that  are  closer  to  each  other  in  position  yield  a 
higher  PDOP,  which  will  degrade  performance. 

Another  concern  is  how  the  antenna’s  gain  pattern  affects  geolocation  capabilities.  At 
lower  frequencies,  there  is  not  much  of  an  issue  because  it  is  easier  to  build  antennas 
which  have  a  good  omnidirectional  gain  pattern.  However,  as  frequency  increases  it 
becomes  more  difficult  to  maintain  a  high  gain  in  a  variety  of  directions  [39].  The  front 
end  of  the  sensor  platform  must  have  a  well  engineered  antenna  array  which  provides  a 
gain  pattern  suitable  for  the  specific  mission. 
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2.5  Performance  Measures 


There  are  many  ways  in  whieh  geoloeation  methods  ean  be  analyzed  and  measured  to 
eompare  their  performanee.  Some  methods  rely  on  statistieal  means  to  measure  the  lower 
varianee  bounds  on  the  solution.  Other  methods  evaluate  the  geometry  of  the  problem 
spaee  for  information  about  how  well  the  geoloeation  methods  might  work.  A  simple 
method  to  gauge  performanee  of  an  estimator  is  to  run  many  iterations  of  a  simulation  and 
traek  the  average  error  over  the  ensemble. 

2.5.1  Cramer  Rao  Lower  Bound. 

The  CRLB  is  a  statistieal  lower  bound  on  the  varianee  of  an  unbiased  estimator  [40]. 
Comparing  results  using  the  CRLB  provides  insight  into  the  effeetiveness  of  a  partieular 
algorithm,  and  ean  aid  in  the  deteetion  of  biases  and  ineflieieneies.  The  CRLB  is  a 
funetion  of  several  faetors,  ehief  among  them  the  number  and  geometry  of  the  sensors,  the 
measurement  type,  and  the  dimension  of  the  solution  spaee  (i.e.  2D  or  3D). 

At  the  heart  of  the  CRLB  is  the  Fisher  Information  Matrix  (FIM).  The  FIM  measures 
the  amount  of  information  that  a  random  variable  X  holds  in  relation  to  an  unknown 
parameter  6,  aeeording  to  the  following  equation  [41]: 


FIMiO)  =  -E 


d'^ln(p(X;e)) 


(2.42) 


The  expeeted  value  of  the  seeond  derivative  of  the  signal’s  probability  distribution 
funetion  (PDF)  provides  the  FIM  for  eaeh  independent  observation.  Another  useful 
feature  of  the  FIM  is  that  eaeh  observation’s  FIM  ean  be  added  together  to  assess  arrays  of 
signal  samples,  so  long  as  the  observations  are  independent. 

In  a  general  ease,  this  researeh  assumes  that  the  distribution  of  the  noise  summed 
with  the  signal  is  Gaussian,  and  so  has  the  following  PDF,  p(X;  6)  [41]: 
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where  cr^  is  the  variance  of  each  observation,  Ns  is  the  number  of  samples  collected,  X  is 
the  received  signal,  and  Sr  is  the  transmitted  signal.  Taking  the  natural  logarithm  then 
differentiating  (2.43)  twice  with  respect  to  the  variable  of  interest  6  results  in  the  following 
equation  [41]: 

_  y  [  dSr[n-fi]'^ 

^2  ^  [  de  J 

— lnp(X;  6)  =  -  (2.44) 

(7-2 

In  order  to  find  the  CRLB,  the  result  of  (2.44)  is  inverted  as  follows  [41]: 
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The  generalized  CRLB  formula  in  (2.44)  and  (2.45)  will  change  based  on  the  relationship 
between  Sr  and  6. 

2.5.2  Dilution  of  Precision  (DOP). 

An  important  factor  to  consider  in  geolocation  accuracy  is  the  geometry  between  the 
satellites  and  the  transmitter.  A  large  intersection  angle  is  desirable  [8],  and  increasing  the 
distance  between  satellites  is  helpful  for  this.  There  is  a  limit  to  how  far  apart  the  sensors 
can  be  positioned  before  performance  is  degraded.  Ho  and  Chan  provide  a  set  of  equations 
that  can  be  used  to  determine  this  limit  [8].  Also,  there  is  a  balance  that  should  be 
considered  between  the  PDOP  and  the  target  transmitter’s  beamwidth.  In  order  to 
maximize  the  SNR,  satellites  must  remain  within  one  beamwidth  of  the  transmitting 
antenna.  The  transmitter’s  beamwidth  is  a  concern  in  designing  the  orbit  used  in 
transmitter  detection.  Some  applications  cannot  assume  a  known  transmitter  beamwidth 
and  transmitter  bearing,  and  would  need  to  keep  the  range  between  the  satellites  small  to 
ensure  that  the  signal  is  collected  by  all  of  the  sensors  simultaneously.  The  scenario 
necessitates  a  poor  dilution  of  precision  in  order  to  meet  the  required  signal  collection 
parameters  for  the  geolocation  algorithms  to  run. 

PDOP  is  a  concern  in  all  localization  methods,  but  can  be  much  more  pronounced  in 
far- field  applications.  The  relative  distance  between  a  sensor  and  the  transmitter  is 
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generally  much  greater  than  the  distance  between  the  sensors,  which  makes  the  various 
range  estimates  from  sensors  to  transmitter  very  close  numerically.  PDOP  is  calculated  by 
a  root-sum-square  of  the  autocorrelation  values  of  distance  measurements  in  the  sensor’s 
relation  to  the  transmitter.  The  PDOP  equation  is  as  follows  [42] : 

PDOP  =  ^crl  +  ^2  +  crj  (2.46) 

where  cr^,  (Ty,  and  cr^  represent  position  uncertainty  in  their  respective  coordinate 
directions.  These  values,  referred  to  as  “sigma  values”,  increase  as  the  range 
measurements  between  sensors  increase  in  correlation.  Larger  sigma  values  yield  more 
diluted  position  estimates  due  to  increased  overlap  in  noisy  range  estimates.  Figure  2.6 
shows  examples  of  both  good  and  poor  PDOP. 


Poor  GDOP 


Well  spaced  satellites-  low 
uncertainty  of  position 


Good  GDOP 


Poorly  spaced  satellites  - 
high  uncertainty  of  position 


Geometry  in  2-D  (GPS 
Basics,  2000} 


Figure  2.6:  Poor  PDOP  can  increase  uncertainty  in  the  positioning  accuracy  of  geolocation 
systems,  while  good  PDOP  can  lessen  the  uncertainty  [43]. 


2.6  Cubesat  Characteristics 

The  Cubesat  standard  is  a  relatively  recent  development  in  the  satellite  community. 
In  the  1990’s,  many  satellite  programs  began  designing  increasingly  large  and  complex 
satellites  which  were  envisioned  to  save  money  by  doing  many  missions  onboard  rather 
than  just  one.  Many  programs,  however,  found  that  the  increased  complexity  of  the 
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platform  actually  led  to  more  problems,  more  delays,  and  more  cost  [32].  More  complex 
missions,  combined  with  the  ever-shrinking  microprocessor,  led  to  a  new  trend  to  simplify 
and  miniaturize  satellites.  Smaller  and  cheaper  was  perceived  as  better,  and  thus  a 
standard  was  created.  The  Cubesat  is  a  satellite  standard  started  by  Puig-Suari  and  Twiggs 
at  Cal  Poly  and  Stanford  in  1999  [44].  The  base  platform  size  is  lOx  lOx  10  cm^  and 
weighs  no  more  than  1.33  kg  -  this  is  called  lU  (One  Unit).  For  the  past  decade,  Cubesats 
have  been  primarily  used  for  university  research  projects  due  to  their  small  payload  and 
affordability.  Educational,  military,  and  commercial  interest  has  grown  substantially  as 
their  demonstrated  capabilities  have  expanded. 

The  advantages  of  using  satellites  to  perform  a  geolocation  mission  stem  from  their 
vantage  point  from  orbit.  The  field  of  view  is  much  greater  than  what  could  be  achieved 
by  an  airborne  platform.  Additionally,  there  are  advantages  in  policy,  as  there  are  no 
overfly  restrictions  for  satellites.  The  likelihood  of  a  satellite  being  forcibly  removed  from 
service  are  low  as  well,  so  a  satellite  is  a  resilient  host  for  the  geolocation  mission. 

The  idea  of  developing  space-based  assets  on  a  short  timeframe  and  a  small  budget 
has  been  a  popular  topic  in  recent  years.  The  ORS  Office  has  the  following  task  [45]: 

...plan  and  prepare  for  the  rapid  development  of  highly  responsive  space 
capabilities  that  enable  delivery  of  timely  warfighting  effects  and,  when 
directed,  develop  and  support  deployment  and  operations  of  these  capabilities 
to  enhance  and  assure  support  to  Joint  Force  Commanders'  and  other  users' 
needs  for  on-demand  space  support,  augmentation,  and  reconstitution. 

Cubesats  may  be  an  excellent  way  to  to  effectively  and  efficiently  provide  reconnaissance 
and  intelligence  to  warfighters  on  the  short  timeline  desired  by  ORS.  The  use  of  Cubesats 
specifically  for  this  type  of  mission  is  beneficial  because  of  their  cost  effectiveness  and 
relatively  short  design  and  build  time.  The  low  cost  of  production  for  a  Cubesat  could 
allow  a  set  of  satellites  to  be  built  and  placed  on  the  shelf  for  short-notice  surveillance 
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requirements.  The  time-to-launch  factor  of  a  solution  like  this  is  then  limited  only  by  the 
time  required  to  determine  an  optimal  orbit  and  RF  front  end,  then  launch.  The  ORS 
timeline  goal  of  months  to  a  year  [45]  is  feasible. 

2.7  Cubesat  Capabilities 

Selva  and  Krejci  [32]  give  a  good  introduction  to  the  capabilities  and  limiting  factors 
of  Cubesats  used  in  an  Earth  observation  mission.  One  such  limiting  factor  is  the 
capability  for  data  transmission,  which  is  much  more  restrictive  than  the  capability  for 
data  storage.  A  typical  Cubesat  has  a  downlink  data  rate  of  0.25  Mbps,  and  only  has  about 
5  minutes  of  access  to  the  ground  station  per  pass.  The  short  access  window  leads  to  a 
maximum  transmission  capability  of  9.3  MB  per  pass,  far  below  the  gigabytes  of  data 
storage  that  may  be  available  [32] .  The  difference  could  affect  the  geolocation  mission  if  a 
large  amount  of  data  is  required  to  be  post-processed  at  the  ground  station.  The  amount  of 
computation  that  can  be  performed  on-board  the  satellite  can  mitigate  the  need  to  transmit 
extensive  data  back  to  the  ground  station,  but  at  the  cost  of  power,  computation  time,  and 
required  on-board  storage. 

Another  limiting  factor  which  has  an  immediate  implication  for  geolocation 
algorithm  accuracy  is  the  attitude  determination  sensors.  Current  technology  has  achieved 
accuracies  of  less  than  2°  on  a  deployed  satellite  (CanX-2)  [32].  This  system  uses  sun 
sensors  and  magnetometers  for  attitude  awareness.  A  star- tracker  could  be  employed  with 
better  results,  but  there  doesn’t  appear  to  be  an  operational  Cubesat  with  such  a  system  in 
place  [32].  There  are  currently  star  trackers  in  development  which  would  achieve  attitude 
measurement  accuracy  of  0.5°  [32].  This  accuracy  would  lower  position  uncertainty  on 
the  ground  to  4.5  km  from  17.5  km  for  a  2°  position  knowledge  accuracy. 

The  dimensions  of  the  Cubesats  also  present  some  limitations.  There  are  obvious 
difficulties  that  come  with  a  small  payload,  to  include  less  computational  ability,  less 
robust  communication,  and  shorter  lifespan.  Receiver  capability  can  be  degraded  as  well. 
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The  maximum  diameter  of  an  RF  antenna  on  a  lU  Cubesat  is  10  cm,  which  corresponds 
to  the  wavelength  of  a  3  GHz  signal.  Any  fixed  antenna  put  on-board  the  lU  Cubesat  will 
not  be  a  full  wavelength  for  signals  below  the  C-band.  The  length  restriction  could  affect 
gain  and  SNR  of  the  antenna,  and  may  reduce  communication  or  sensing  capabilities. 
Deployable  antennae  can  be  included  in  the  payload,  but  this  increases  the  risk  of 
mechanical  failure  and  impacts  the  space  available  for  other  hardware  on-board  the 
Cubesat.  This  research  focuses  on  a  6U  Cubesat  architecture,  which  will  lessen  the  effect 
of  such  limitations.  A  6U  Cubesat  has  a  maximum  of  30  cm  length  for  a  fixed  antenna, 
allowing  communications  in  the  L-band. 

2.8  Geolocation  System  Designs 

The  state  of  the  art  in  geolocation  systems  is  evolving  to  include  satellites  as  a  sensor 
platform.  There  are  a  number  of  satellite  systems  which  have  been  designed  to  perform  a 
geolocation  mission.  The  Satelllite  Mission  for  Swarming  and  Geolocation  (SAMSON) 
project  is  a  three-Cubesat  constellation  mission  concept  designed  to  perform  autonomous 
cluster  flight  and  geolocation  of  a  cooperative  RF  emitter  on  the  surface  of  Earth  [46].  The 
SAMSON  project  is  led  by  the  Technion  (the  Isreali  Institute  of  Technology)  and 
supported  by  the  Israeli  space  industry.  As  of  this  writing,  SAMSON  is  in  the  design 
phase,  with  anticipated  launch  in  2015.  The  satellite  swarm  will  use  TDOA  and  FDOA  to 
perform  a  geolocation  estimate  for  use  in  search  and  rescue  operations.  The  system  is 
intended  for  a  ‘lost  at  sea’  scenario  where  a  transmitting  wristband  will  be  worn  by 
boaters  to  provide  the  satellites  a  priori  information  of  the  signal  structure.  The 
geolocation  system  is  expected  to  determine  TDOA  measurements  within  100 
nanoseconds  and  FDOA  measurements  within  0.5  Hertz.  Knowledge  of  satellite 
orientation  and  position  is  expected  to  be  within  0.5  degrees  of  truth  using  a  combination 
of  Earth  and  Sun  sensors,  magnetometers,  and  magnetorquers. 
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Although  Cubesats  are  beginning  to  see  more  interest  for  geolocation  missions, 
larger  satellites  are  predominantly  used.  A  recent  patent  from  Northrop  Grumman 
Systems  Corporation  shows  a  design  for  an  interferometer  based  geolocation  system.  The 
satellite  would  carry  an  interferometer  array  to  produce  a  set  of  direction  of  arrival  vectors 
[47].  The  system  would  only  require  a  single  satellite  to  perform  the  methods  proposed  in 
the  patent.  Qualcomm  Incorporated  have  proposed  a  similar  application,  where  two 
satellites  are  used  to  determine  the  position  of  a  user.  The  patent  filed  by  Levanon  et  al. 
uses  TDOA  and  FDOA  to  perform  geolocation  to  locate  a  cooperative  signal  source. 

2.9  Summary 

This  section  presented  a  background  of  the  geolocation  problem  and  its  solution 
methods.  There  are  two  primary  classical  methods  which  have  been  developed  to  support 
the  research  into  localization  methods.  TDOA  uses  the  difference  in  reception  times 
between  sensors  to  estimate  the  source  of  an  RF  transmission.  AOA  methods  calculate 
lines  of  bearing  from  the  receivers  to  create  a  position  estimate  based  on  where  those 
LOBs  cross.  Another  classical  method,  FDOA,  is  often  used  to  supplement  TDOA  or 
AOA  with  frequency  difference  measurements  if  the  receivers  and/or  the  transmitters  are 
moving.  DPD  is  a  new  localization  method  which  uses  a  cost  function  for  a  single  step 
solution.  DPD  has  been  shown  to  significantly  improve  geolocation  performance  in 
systems  with  a  low  SNR.  The  research  with  respect  to  the  DPD  method  has  been  focused 
on  a  geolocation  problem  measured  in  tens  of  kilometers.  Research  into  the  performance 
of  the  DPD  method  at  standoff  distances  of  hundreds  or  thousands  of  kilometers  is  unique 
to  this  research,  and  is  explored  in  this  thesis. 

Cubesats  are  a  relatively  new  platform  in  the  satellite  community.  Small  in  size  and 
budget,  they  are  ideal  for  programs  where  the  mission  is  computationally  inexpensive  and 
design  time  is  short.  In  a  scenario  where  information  about  a  target  C3  node  is  needed  on 
short  notice,  Cubesats  can  be  deployed  quickly  if  there  is  a  model  in  place  to  determine 
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the  best  solution  for  a  particular  set  of  parameters.  A  pre-built  Cubesat  constellation  could 
be  uploaded  with  the  appropriate  algorithms  in  software  and  deployed  with  little 
preparation  time. 


42 


III.  Methodology 


This  chapter  provides  the  reader  with  a  detailed  description  of  the  simulations 

conducted  and  the  algorithms  used  in  this  work.  The  simulation  environment  is 
described,  to  include  the  transmitter  characteristics,  signal  of  interest  parameters,  and 
methods  used  to  geolocate.  The  performance  metrics  used  to  evaluate  the  various 
geolocation  methods  are  detailed  as  well. 

3.1  Problem  Overview 

The  problem  considered  is  the  geolocation  of  a  fixed  RF  transmitter  with  a  variable 
number  of  Cubesats  in  LEO  above  the  target.  The  satellites  are  capable  of  receiving  the 
signal  while  in  orbit,  and  are  able  to  communicate  with  other  satellites  in  the  constellation 
to  share  the  gathered  information  about  the  received  signal.  The  satellites  are  in  radio 
communication  with  a  ground  station  or  a  network  of  ground  stations  to  communicate 
location  estimates  and  error  parameters  to  the  user.  Figure  3.1  depicts  an  anticipated 
mission  scenario. 

3.2  Simulating  the  Space  Environment 

In  order  to  simulate  the  space  environment  in  the  scenario,  care  has  be  taken  to 
ensure  the  space-borne  sensors  are  modeled  as  moving  according  to  realistic  orbit 
parameters.  Also,  the  simulation  of  the  signal  must  be  modeled  realistically  in  its 
propagation  through  free-space  from  the  transmitter’s  position  to  the  satellites’  positions. 
The  satellites  are  modeled  to  be  in  constant  motion  along  their  orbits,  so  a  Doppler  shift  is 
simulated  and  associated  with  the  signal  as  well. 

The  AGI  System  Tool  Kit  (STK)  is  used  in  this  research  to  simulate  the  satellites’ 
positions  throughout  their  orbits.  STK  is  a  tool  commonly  used  by  astronautical  engineers 
and  satellite  designers  to  build  realistic  scenarios  involving  spacecraft,  heavenly  bodies. 


43 


Unknown  RF 
source 


Communicate 

between 

satellites 


r  \<>r\ 


Uplink  ' 
commands 


Downlink 
L  payload  data 
\  andSOH 


Network  To  and 


Figure  3.1:  Scenario  Overview 


and  Earth-based  objects.  The  tool  allows  for  simulation  of  communication  connections 
between  satellites  and  sensors  on  the  ground,  and  is  able  to  accurately  determine  the 
position  of  each  satellite  when  the  communication  link  is  available.  This  data  can  be 
exported  from  STK  and  then  be  imported  into  Matlab  for  use  in  the  geolocation 
simulation. 

In  a  real  world  scenario,  the  orbit  parameters  would  be  chosen  to  fit  the  mission  -  i.e. 
to  maximize  time  over  target,  number  of  passes  per  day,  etc.  For  this  research,  the  site  of 
the  simulated  transmitter  is  chosen  to  be  AFIT,  at  the  geodetic  coordinates  [39.782272°  N 
-84.082893°  W  0].  The  altitude  of  0  meters  indicates  that  the  transmitter  is  on  the  surface 
of  the  Earth. 

Table  3.1  shows  the  orbit  parameters  for  a  four  satellite  scenario.  When  multiple 
satellites  are  simulated,  the  inclination  angle  is  changed  by  0.5°,  and  the  altitude  is 
increased  by  1  km  for  each  additional  satellite  to  reduce  the  risk  of  collision  in  orbit.  The 
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Table  3.1:  Orbit  Parameters 


Cubesat  1 

Cubesat  2 

Cubesat  3 

Cubesat  4 

Altitude 

450  km 

451  km 

452  km 

453  km 

Eeeentrieity 

0“ 

0° 

0° 

0° 

Inelination 

o 

o 

49.5° 

49° 

50.5° 

Period 

93.6  minutes 

93.6  minutes 

93.6  minutes 

93.6  minutes 

Maximum  Pass  Time 

10  minutes 

10  minutes 

10  minutes 

10  minutes 

details  of  the  orbit  parameter  trade  spaee  are  not  the  foeus  of  this  researeh,  but  they  are 
ehosen  to  inerease  time-over-target  and  improve  eonstellation  geometry  while  over  the 
target.  While  over  the  target,  the  geometry  ehanges  from  a  point  with  one  satellite  to  a  line 
with  two  satellites,  a  triangle  shape  with  three  satellites,  or  a  parallelogram  shape  with 
four.  The  PDOP  for  these  formations  is  better  than  that  offered  by  a  straight  line 
formation.  Figure  3.2  illustrates  the  seenario  from  Table  3.1,  with  four  satellites  orbiting 
over  the  target  in  three  dimensions  in  STK. 


Figure  3.2:  3D  Orbit  Over  Target 
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The  simulation  considered  in  this  research  only  uses  a  single  pass  over  the  target  to 
generate  the  geolocation  estimate.  The  Cubesats  receive  the  signal  for  515  seconds  along 
the  orbit,  and  each  second  the  satellite  positions  are  updated  according  to  the  orbital 
information  generated  by  STK.  Figure  3.3  shows  the  considered  orbit  over  the  target.  The 
pass  used  in  this  research  does  not  pass  directly  over  the  target,  nor  is  it  the  longest  pass 
possible.  In  realistic  scenarios,  it  would  be  unlikely  to  pass  directly  over  the  target,  as  its 
position  is  unknown.  The  closest  point  in  the  satellite’s  trajectory  is  approximately  770 
km  from  the  target. 


Figure  3.3:  The  portion  of  the  orbit  considered  in  this  research  does  not  pass  directly  over 
the  target.  The  red  star  represents  the  target  and  the  colored  asterisks  represent  the  Cubesat 
positions  over  time. 


3.3  System  Parameters 

The  geolocation  system  is  comprised  of  N  satellites.  The  composition  of  the 
constellation  is  varied.  The  satellite  is  modeled  as  having  four  antennas  with  associated 
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receiver  radios  on-board  each  of  the  satellites,  and  the  radios  are  able  to  receive  the  SOI 
and  accurately  detect  frequency  and  phase  information.  Reference  timing  on  a  Cubesat  is 
typically  based  on  GPS  because  there  is  not  enough  space  for  a  finely  calibrated  oscillator 
on-board.  There  are  some  small  oscillators  available,  such  as  the  Symmetricom  Quantum 
[48]  and  the  Jackson  Labs  CSAC  [49].  The  difficulty  with  including  these  oscillators  is 
that  they  would  necessitate  synchronization  between  sensors,  which  adds  further  error 
potential  and  requires  more  power.  For  this  research,  therefore,  it  is  assumed  that  GPS 
timing  is  used,  and  that  the  timing  has  a  Gaussian  error  with  mean  zero  and  a  standard 
deviation  of  40  nanoseconds  [35]. 

Once  the  signal  is  received,  sampled,  and  filtered,  it  can  be  sent  to  the  on-board 
processor.  This  research  assumes  that  the  geolocation  algorithms  are  executed  on-board 
the  satellite.  With  the  sampled  signal  available,  the  processor  can  extract  the  information 
required  to  perform  geolocation.  As  discussed  in  Section  2.1,  this  may  include 
information  about  timing,  phase,  frequency,  or  some  combination  of  all  of  these.  The 
satellite  then  transfers  the  extracted  data  to  a  central  node.  From  there,  data  are  compiled 
and  an  estimate  is  computed.  The  estimate  may  be  output  to  the  ground  station,  saved  in 
memory  for  tracking  purposes,  or  both.  In  addition  to  the  position  estimate,  an  estimate  of 
the  solution’s  variance  can  be  determined.  The  combination  of  the  estimate  and  its 
variance  provides  the  user  with  insight  into  where  to  search  for  the  target.  Figure  3.4 
shows  an  overview  of  the  processing  flow  from  signal  reception  to  final  position  estimate. 

3.4  Models  for  Signal  Generation,  Transmission,  and  Reception 

The  signal  produced  by  the  transmitter  is  assumed  to  be  a  narrowband  signal, 
transmitted  with  enough  power  to  remain  above  the  noise  floor  when  received  at  LEO. 

The  bandwidth  of  the  narrowband  signal  is  approximately  100  kHz.  The  waveform 
generated  in  the  simulation  is  a  50  microsecond  square  pulse,  as  shown  in  Figure  3.5. 
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Figure  3.4:  Processing  Flow 
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Figure  3.5:  Simulated  Transmission  Signal 


The  signal  can  be  broadcast  with  varying  power,  which  changes  the  quality  of  the 


received  signal.  The  simulated  transmission  signal  is  comprised  of  500  sinusoids  added 
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together  to  form  a  square  pulse.  The  equation  used  in  the  simulation  to  model  the 
transmission  signal  s(t)  is  the  following: 

500 

s(t)  = '^a„cos  (Innfct)  (3.1) 

n=l 

where  is  the  amplitude,  fc  is  the  carrier  frequency  of  the  signal,  and  t  is  the  time  vector. 


For  a  square  pulse,  the  amplitude  is  a  sine  function,  which  is  modeled  as  follows: 


Qn 


2Asin(2;7m|j 

2nn 


(3.2) 


where  A  is  the  amplitude  of  the  pulse,  k  is  the  width  of  the  pulse,  and  T  is  the  period.  To 
model  channel  effects  as  the  signal  is  broadcast,  s{t)  is  delayed  by  t,  shifted  by  the 
Doppler  frequency  fd,  and  attenuated  by  a,  yielding  the  signal  r{t).  The  received  signal  is 
modeled  as  follows: 


r{t)  =  a 


500 


^  Qncos  (2nnfc  (t  -  r)) 


L  n=\ 


oA^fdt 


+  w{t) 


(3.3) 


where  a  noise  term  w{t)  is  included  in  r{t)  to  simulate  additive  white  Gaussian  noise 
(AWGN)  that  would  be  present  in  the  received  signal.  The  received  signal  can  vary  widely 
depending  on  the  variables  inserted  into  (3.3).  The  noise  added  to  the  signal  is  of  a 
magnitude  dictated  by  the  SNR  at  the  receiver.  A  signal  with  a  large  SNR  may  be  from  a 
powerful  radar  system,  while  a  signal  with  a  small  SNR  may  be  attributed  to  a  smaller 
Earth-based  communication  system. 

The  time  shift  r  that  is  applied  to  the  signal  corresponds  with  the  time  delay 
associated  with  LOS  signal  propagation.  It  is  assumed  in  this  research  that  the  signal 
propagates  at  c  =  2.9978  x  10*  m/s.  A  signal  traveling  directly  from  the  transmitter  to  the 
Cubesat  directly  overhead  at  an  altitude  of  450  km  is  assumed  to  have  a  time  delay  r  of 
1.5  milliseconds. 

The  Doppler  shift  corresponds  with  the  velocity  of  the  satellite  relative  to  the 
transmitter.  The  typical  speed  of  a  satellite  in  LEO  is  on  the  order  of  7,500  m/s.  This 
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yields  a  substantial  Doppler  shift  which  can  be  measured  to  provide  information  about  the 
emitter’s  position.  The  speed  is  related  to  the  position  of  the  receiver  as  well  as  the  timing 
of  the  samples.  The  Cubesat  orbit  position  error  is  generally  less  than  10  meters  when 
using  an  on-board  star  tracker  [50].  A  Gaussian  error  with  mean  zero  and  a  standard 
deviation  of  10  meters  is  added  to  each  receiver  position  to  simulate  these  errors.  Using 
the  equation  for  Doppler  shift  for  a  stationary  transmitter  and  a  moving  observer,  the 
Doppler  shifted  frequency  fd  is  as  follows: 


where  Vr  is  the  satellite’s  velocity  relative  to  the  transmitter  and  fc  is  the  carrier  frequency 
of  the  transmitted  signal.  If  the  satellite  is  moving  toward  the  transmitter,  the  relative 
velocity  will  be  positive  and  fd  will  be  larger  than  If  the  satellite  is  moving  away  from 
the  transmitter,  the  relative  velocity  will  be  negative  and  fd  will  be  smaller  than  fc. 

It  should  be  noted  that  atmospheric  conditions  may  change  some  of  the  signal 
propagation  parameters,  such  as  the  speed  and  path  loss  of  the  propagated  signal. 

Applying  a  more  robust  model  for  atmospheric  conditions  will  improve  the  accuracy  of 
the  simulations,  but  accounting  for  atmospheric  attenuation  characteristics  is  not  in  the 
scope  of  this  research.  Accounting  for  these  characteristics  is  an  area  for  further  research. 

3.5  Algorithms 

Four  algorithms  are  considered  in  this  research  for  the  estimation  of  transmitter 
position.  Two  classical  methods  -  TDOA  and  AOA  -  are  implemented  along  with  a  version 
of  the  DPD  method  and  the  IRF  method.  These  algorithms  were  developed  in  Chapter  II. 
The  following  sections  will  detail  their  use  in  the  simulations  built  for  this  research. 

3.5.1  Time  Difference  of  Arrival. 

Chapter  II  describes  the  overall  goal  of  the  TDOA  method  -  to  extract  meaningful 
timing  information  from  a  signal  and  use  that  information  to  obtain  the  difference  in  time 
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of  arrival  between  a  group  of  sensors.  Once  the  TDOA  information  is  extracted,  a 
multilateration  algorithm  follows  to  obtain  a  position  estimate.  A  closed  form 
multilateration  approach  formulated  by  Ho  and  Chan  [8]  is  used  in  this  research.  The 
algorithm  is  developed  as  presented  by  Kulumani  [29],  and  as  described  in  Sections 
2.3.3. 1  to  2. 3. 3. 2.  The  closed  form  approach  is  chosen  because  it  will  always  converge, 
though  it  may  be  more  computationally  burdensome. 

An  important  change  to  Kulumani’s  method  must  be  implemented  in  this  simulation. 
Kulumani  based  his  calculations  on  a  spherical  Earth,  and  this  simulation  uses  an  elliptical 
Earth  model,  the  WGS  coordinate  system,  to  determine  the  positions  of  the  transmitter 
and  sensors.  The  radius  of  the  Earth  which  is  used  in  (2.26)  must  be  modified,  or  else  the 
resulting  emitter  position  will  have  a  significant  source  of  error  built  into  the  equation. 
This  error  will  increase  as  the  transmitter  is  moved  away  from  the  equator,  so  the  error  is 
maximized  when  the  transmitter  is  located  at  one  of  the  Earth’s  poles.  It  is  prudent, 
therefore,  to  determine  the  approximate  Earth  radius  at  the  latitude  of  the  transmitter.  In 
practice,  the  latitude  may  not  be  implicitly  known.  If  a  general  target  area  is  known  (such 
as  if  it  is  known  that  the  transmitter  is  located  near  a  certain  city),  the  latitude  could  be 
approximated  by  the  user.  If  the  target  whereabouts  are  completely  unknown,  the  process 
could  be  conducted  iteratively,  where  Earth’s  equatorial  radius  is  used  first  and  then  the 
latitude  is  refined  after  the  estimate  is  computed.  This  research  assumes  that  the 
transmitter  is  known  to  be  in  or  near  a  large  city,  whose  approximate  latitude  is  known. 
The  equation  to  find  the  radius  of  the  Earth  re(i/r)  at  a  given  geodetic  latitude  is  the 
following  [51]: 


a(iA)  = 


(a^  cos  +  {b^  sin  ipY 


(3.5) 


(a  cos  iff)^  +  (b  sin 

where  ij/  is  the  geodetic  latitude  of  the  transmitter,  a  is  the  Earth’s  equatorial  radius  (6378 
km),  and  b  is  the  Earth’s  polar  radius  (6356  km).  Using  the  calculated  value  of  from  3.5 
should  reduce  a  significant  error  source,  relative  to  that  which  would  be  found  when  using 
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Kulumani’s  method.  The  TDOA  method  can  only  be  used  when  three  or  more  satellites 
are  available  in  the  geolocation  system. 

3.5.2  Angle  of  Arrival. 

The  angle  of  arrival  is  calculated  using  the  MUSIC  algorithm,  which  was  developed 
in  Section  2.3. 1.4.  Simulations  associated  with  this  research  assume  a  single  signal  is 
present,  even  though  the  MUSIC  algorithm  is  applicable  in  multi-signal  environments 
(provided  there  are  more  antennas  than  the  number  of  signals  present).  The  satellite  is 
simulated  to  have  an  antenna  array  consisting  of  two  antennas  per  direction  estimated. 
Azimuth  and  elevation  angles  are  required  to  construct  a  three  dimensional  line  of 
bearing.  Once  the  bearing  angles  are  estimated,  the  line  of  bearing  is  calculated  beginning 
at  the  satellite’s  coordinates  and  propagating  out  in  the  direction  of  the  emitter.  LOB’s  are 
defined  using  a  point  of  origin  (i.e.  the  satellite’s  position)  and  an  associated  bearing 
angle.  As  more  LOB’s  are  created,  they  can  be  aggregated  to  find  the  origin  of  the  signal 
by  determining  cross  points.  There  will  be  multiple  points  where  the  lines  cross  due  to 
noise  in  the  signal.  When  three  LOB’s  are  present  in  2D  this  will  create  a  ‘cocked  hat’,  as 
described  by  Stein  [31]. 

This  research  presumes  a  3D  coordinate  system,  which  makes  the  process  of  finding 
crossings  in  LOB’s  more  complicated.  A  LOB  in  3D  has  an  origin  point  and  two 
associated  angles  for  azimuth  and  elevation,  (6,  f).  Multiple  LOB’s  from  noisy  signals 
may  never  actually  intersect  in  three  dimensions.  Instead  of  finding  the  center  of  Stein’s 
‘cocked  hat’,  a  better  approach  is  to  use  a  least  squares  method.  Least  squares  essentially 
finds  a  position  estimate  in  3D  space  which  minimizes  the  distance  between  that  point  and 
each  of  the  3D  LOBs.  Least  squares  is  used  to  obtain  the  final  estimate  of  the  target’s 
position  through  the  AOA  method,  as  described  in  Section  2. 3. 2.2. 

When  three  or  more  sensors  are  available  to  geolocate  in  the  simulation,  a  3D 
estimate  is  created  for  each  time  epoch  when  a  measurement  is  made.  When  only  one  or 
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two  satellites  are  present,  the  LOBs  are  aggregated  over  a  period  of  time,  and  a  least 
squares  estimate  is  found  using  (2.13). 

3.5.3  Position  Determination  from  LOB’s. 

Once  the  LOB’s  are  calculated  by  the  MUSIC  algorithm  described  in  Section  2.3. 1.4, 
they  are  combined  with  the  associated  satellite  positions  to  determine  the  target  estimate 
through  least  squares  approximation.  The  bearing  lines  are  defined  by  two  points  in  space 
-  an  origin  and  an  end  point.  The  origin  is  the  satellite’s  position,  v,  =  [jc;  y,  nY ,  and  the 
end  point  w,  =  ^Xprop  yprop  Zprop\  is  found  by  propagating  the  angle  of  arrival  out  from  the 
satellite’s  position.  The  range  of  propagation  should  be  far  enough  to  ensure  the  surface  of 
the  Earth  is  surpassed.  For  this  simulation  the  propagation  range  Vprop  is  calculated  as 
follows: 

rprop  =  1.5  X  argmax(r;)  (3.6) 

N 

where  r,  is  the  Euclidean  range  from  the  satellite  to  the  target. 

3.5.4  Direct  Position  Determination. 

The  DPD  method  is  based  upon  classical  TDOA  and  FDOA  methods,  but  formulates 
emitter  positions  in  a  single  step  rather  than  the  dual  steps  of  calculating  timing  and 
Doppler  prior  to  estimating  position.  The  original  method  proposed  by  Weiss  [18] 
requires  the  computations  to  be  done  at  a  single  node,  and  requires  a  series  of  CAE’s  to  be 
computed.  Pourhomayoun  and  Fowler  have  simplified  the  computations  and  distributed 
the  processing  in  their  paper  [16].  Weiss’  original  method  is  used  in  this  research. 

The  signal  model  for  Weiss’  method  is  more  complex  than  the  signal  models  used  for 
the  TDOA  discussion  in  Section  2. 3. 1.1.  The  received  signal  at  the  sensor  s„[zn]  is 
defined  as  follows  [18]: 

s„[zw]  =  h„a„(p)5,„[zw]  +  w„[zw],  zw  6  [0, . . . ,  A,  -  1}  (3.7) 

where  s™  is  a  1  x  vector,  is  a  complex  number  representing  attenuation,  a„(p)  is  the 
array  response  at  the  sensor  from  the  transmitter’s  coordinates  p,  Stn  is  the  transmitted 
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signal,  w„[m]  is  AWGN,  m  is  the  discrete  time  sample  number,  and  Ng  is  the  number  of 
samples.  In  order  to  extract  the  useful  information  from  the  received  signal  waveform,  a 
discrete  Fourier  transform  (DFT)  is  performed  on  the  received  signal.  The  frequency 
domain  signal  s„  follows  [18]: 

s„[m]  =  +  w„[m],  m  e  (0, . . . ,  A,  -  1}  (3.8) 

where  s„  is  a  1  x  vector,  <x>„  =  ^,k\s  between  0  and  -  1,  r„(p)  is  the  propagation 
delay,  and  to  is  the  transmit  time.  To  optimize  the  process,  the  equation  for  the  cost 
function  2(p)  is  minimized  as  follows  [18]: 

N  N  Ns-l  ^ 

Q(P)  =  Yj  IIs™I|'  -  ^  af(p)  £  e^'"‘[^"(P)"^«]4[m]s™[m]  (3.9) 

n=l  n=l  k=Q 

The  equation  for  2(p)  is  cumbersome,  but  it  can  be  manipulated  by  searching  for  the 
maximum  of  the  second  term  rather  than  the  minimum  of  the  entire  expression.  Searching 
for  the  maximum  may  delete  a  term  in  the  equation,  but  there  is  still  an  unknown, 
immeasurable  quantity  in  the  equation,  s*^.  Working  under  the  assumption  that  the 
original  signal’s  waveform  and  transmit  time  are  unknown,  this  term  must  be  resolved.  In 
order  to  accomplish  this,  the  vectors  d„  and  s  are  defined  as  follows  [18]: 

dn  =  [dniO),...,dniN,-l)f  (3.10a) 

4(k)  =  (3.10b) 

s  =  . . . ,  StniNs  -  ]  (3. 10c) 

The  cost  function  in  (3.9)  can  then  be  rewritten  as  follows  [18]: 

N 

e(p)  =  £isXP  (3.11) 

n=l 

The  transmitted  signal  s  is  independent  of  the  sensors,  so  if  a  matrix  D  is  defined  as 

N 

D  =  Z  dfdn,  2(P)  can  be  expressed  as  follows  [18]: 

n=i 

e(p)  =  s^Ds  (3.12) 
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The  vector  s  is  unknown,  but  it  can  be  assigned  to  be  the  eigenvector  corresponding  to  the 
largest  eigenvalue  in  D.  This  assignment  achieves  the  maximum  of  2(p).  The  only 
unknown  in  (3.12)  is  D,  because  s  is  a  fixed  and  unknown  quantity.  Therefore,  2(p)  is 
maximized  by  defining  the  vector  s  as  the  eigenvector  corresponding  to  the  largest 
eigenvalue  of  the  matrix  D,  an  x  matrix  [18]. 

The  DPD  method  relies  on  a  grid  search  model  to  provide  an  estimate  for  the 
emitter’s  position  [18].  The  literature  on  the  DPD  method  primarily  uses  a  two 
dimensional  model,  leading  to  a  square  grid  of  possible  emitter  positions.  Using  satellites 
requires  a  three  dimensional  grid  to  ascertain  the  best  estimate.  The  cube  begins  with  125 
points  spaced  evenly  throughout  the  structure,  then  adds  a  random  error  factor  with  zero 
mean  and  a  standard  deviation  which  starts  at  100  meters  and  shrinks  to  10  meters  as  the 
cube  size  is  decreased.  The  size  of  the  bounding  box  is  200  km  on  each  side  to  start,  and  is 
decreased  by  20  percent  after  each  iteration  until  the  box  is  10  meters  on  each  side. 

3.5.5  Instantaneous  Received  Frequency  Method. 

The  IRF  method  uses  the  changing  velocity  of  the  satellite  relative  to  the  emitter  to 
approximate  the  position  of  the  transmitter.  As  discussed  in  Section  2.3.5,  a  grid  of 
possible  emitter  locations  is  defined,  and  estimated  instantaneous  received  frequency  is 
calculated  for  each  grid  point. 

For  this  research,  the  simulation  creates  a  three  dimensional  cube  with  the  center 
point  at  the  true  emitter  position,  plus  an  error  factor  which  is  random  and  less  than  100 
meters.  The  grid  initially  consists  of  125  points,  or  five  grids  which  are  each  5x5,  spaced 
evenly  throughout  the  cube.  The  size  of  the  bounding  box  is  200  km  on  each  side  to  start. 
Once  the  estimated  IRF  is  calculated  for  each  point  in  the  cube,  the  closest  match  is 
decided  and  the  next  cube  is  constructed  with  that  point  as  the  center.  The  dimension  of 
the  bounding  cube  is  reduced  by  20  percent  until  the  box  is  10  meters  on  each  side. 
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In  order  to  calculate  the  instantaneous  frequency  of  the  received  signal,  the  Welch 
method  is  used  [52].  The  simulation  for  the  IRF  method  uses  a  signal  of  length 
Ns  =  10, 000  samples,  a  Hanning  window  with  a  length  of  1000  samples,  and  an  overlap 
of  500  samples. 

3.6  Performance  Measurements 

The  performance  metrics  used  to  evaluate  the  methods  examined  in  this  simulation 
scenario  are  the  Cramer- Rao  lower  bound  (CRLB)  and  the  position  dilution  of  precision 
(PDOP).  These  metrics  are  related,  but  provide  the  designer  with  different  information 
about  the  target  solution  that  may  be  useful  in  design  decisions.  The  generalized  CRLB  is 
developed  in  Section  2.5.1,  and  is  further  detailed  below  according  to  the  various 
measurement  types  which  are  used  in  this  research. 

3.6.1  Cramer-Rao  Lower  Bound. 

The  CRLB  is  applied  in  this  research  to  show  how  the  various  methods  of 
geolocation  used  in  the  simulation  perform  relative  to  the  statistical  lower  bound  on  the 
variance  of  the  position  estimate.  The  CRLB  is  the  inverse  of  the  FIM  for  the  collected 
data,  and  a  convenient  form  of  the  FIM  for  the  data  f(x)  is  as  follows  [53]: 

J  =  H^C  H  (3.13) 


where  H  is  the  Jacobian  and  C  is  the  covariance  matrix  for  the  collected  signal.  The 
Jacobian  is  calculated  according  to  the  positions  of  the  sensors  in  addition  to  the 
parameters  of  interest.  The  basic  form  of  the  Jacobian  H  can  be  represented  as  the 
following  [53]: 
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where  5i  -  54  are  the  equations  for  the  parameters  of  interest,  and  p  is  the  transmitter 
position  and  frequency,  p  =  [xt  yt  Zt  fc\  ■  The  parameters  which  are  of  interest  in  this 
research  are  timing,  frequency,  and  bearing,  respectively.  Bearing  encompasses  both 
azimuth  and  elevation  angles.  The  equations  for  -  54  are  the  following  [53],  [54]: 


Ri-Rj 


(3.15) 


_  fc  dRi 
^  c  dt 
180 

■^3  =  o:i(az)  + - arctan 

7T 


■^4  ~  ^i{el) 


arccos 


-  yi 

X,  -  Xi 


Zt  -  Zi 


(3.16) 


(3.17) 


(3.18) 


where  Ri  is  the  range  to  from  the  satellite  to  the  transmitter  and  or,-  is  the  direction  of  the 
antenna  main  lobe.  The  number  of  parameters  which  are  used  to  calculate  the  CRLB  for 
position  is  dependent  on  the  geolocation  method  being  used.  For  example,  when  just  the 
DPD  method  is  used,  only  5i  and  $2  are  used  for  timing  and  frequency  information.  The 


elements  of  the  Jacobian  in  (3.14)  can  be  computed  with  (3.15)  -  (3.18)  as  follows: 
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The  matrix  C  has  the  eovarianees  for  the  parameters  of  interest  along  the  diagonal, 


aeeording  to  the  following: 

I  Cl  0  0 

0  •••  0  (3.19) 

0  0  Civ 

where  N  is  the  number  of  parameters  being  used  to  find  the  CRLB,  and  the  parameters  of 
interest  are  assumed  to  be  uneorrelated.  The  diagonal  elements  of  C  ean  be  ealeulated  by 
finding  the  expeeted  value  of  the  parameter  equations  aeeording  to  (3.15)  -  (3.18)  as 
follows: 


^  (var(R;)  +  var(Ry)) 

(3.20) 

(3.21) 

180  \ 

E[s3Sj]  =  + - varie^J 

TT  ^  ' 

(3.22) 

^  180  \ 

E{sas^]=  ai(^ei)+  var(0e/) 

TT  ^  ' 

(3.23) 

Onee  the  Jaeobian  and  the  eovarianee  matrix  are  ealeulated  they  ean  be  used  in  the 
following  equation  to  ealeulate  the  CRLB  for  position  estimation: 


CcRLBiv)  =  r '  =  [h^c-'h] 


(3.24) 


The  following  seetions  will  deseribe  the  varianees  for  the  individual  parameters  of 
interest. 

3.6.1. 1  Time  Delay  / Range  Measurements  CRLB. 

To  measure  the  TDOA  position  estimation  performanee,  the  CRLB  for  time  delay  is 
ealeulated,  then  the  parameters  are  ehanged  to  measure  the  range  estimator  performanee. 
The  CRLB  var(to)  for  time  delay  is  the  following  [41]: 


var(fo)  > 


1 


(3.25) 
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where  Brms  is  the  RMS  bandwidth  of  the  received  signal.  Brms  is  defined  as  follows  [41]: 


B 


rms  — 


f  (2;r/)2|S(/)pJ/ 


A 


f  \s(f)m 


(3.26) 


In  order  to  change  the  CRLB  from  being  defined  with  respect  to  r  to  being  defined  with 
respect  to  range,  a  change  of  variables  is  performed  such  that  R  =  ^,  and 
CRLB^  =  (^)  CRLBfg.  After  this  conversion,  the  final  CRLB  var(R)  for  the  range 
equation  is  as  follows  [41]: 


var(R)  > 


cV4 

SNRxBj^, 


(3.27) 


3.6.1.2  Bearing  Measurements  CRLB. 

Assuming  that  the  wavefront  is  a  plane,  the  scenario  is  illustrated  as  shown  in  Figure 
2.5.  The  sampled  signal  X[n]  at  the  sensor  can  be  written  as  follows  [53]: 


X[n]  =  A  cos  +  w\n\ 


(3.28) 


where  A  is  the  signal’s  amplitude,  Q.s  is  the  spatial  frequency  (rad/sensor),  and  ^  is  a  phase 
delay  factor.  The  parameters  in  (3.28)  which  will  fold  into  the  equation  for  CRLB  are 
6  =  [A  0]^.  The  parameter  of  interest  in  this  case  is  the  AOA  (3,  which  is  a  function  of 
Cls-  The  other  two  parameters  are  considered  nuisance  parameters,  but  must  remain  in  the 
equations  even  though  they  are  not  the  parameter  of  interest.  (3  can  be  written  as  follows 
[53]: 


(3.29) 

where  d  is  the  distance  between  antennas  in  the  receiving  array.  The  Fisher  Information 


Matrix  1(0)  for  this  problem  is  the  following  [41]: 
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(3.30) 
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The  variance  bound  results  from  the  inverse  of  the  Fisher  Information  Matrix  along  with 
the  change  of  variables  in  (3.29).  The  CRLB  vdsifi)  for  a  direction  finding  algorithm  is 
bounded  as  follows  [41]: 


var(^)  > - — ^ -  (3.31) 

(27r)25A^RxM|^(f) 

where  M  is  the  number  of  antennas  in  the  receiving  antenna  array,  L  is  the  total  length  of 
the  array,  and  A  is  the  wavelength  of  the  signal  of  interest. 

3.6.1.3  Frequency  Estimation  CRLB. 

The  lower  bound  on  the  frequency  estimation  for  a  collected  signal  can  also  be  found 
from  the  inverse  of  (3.30),  without  a  change  of  variable.  The  lower  bound  for  frequency 
estimation  is  the  following  [41]: 


var  /o  > 


12 


{lnfSNRxN{m-l) 


(3.32) 


The  frequency  accuracy  decreases  as  so  the  length  of  the  collected  signal  has  an 
important  impact  on  the  performance  of  the  estimator. 

3.6.2  Dilution  of  Precision. 

Dilution  of  precision  measures  the  quality  of  the  solution  with  respect  to  the  position 
of  the  sensors  and  the  transmitter,  as  discussed  in  Section  2.5.2.  This  research  uses  PDOP 
instead  of  GDOP  because  timing  error  found  in  the  GDOP  calculations  is  not  considered 
for  the  algorithms  that  are  investigated  in  this  research.  The  DOP  is  related  to  the  CRLB 
[40].  If  the  estimator  is  efficient,  then  the  relationship  between  the  DOP  and  the  CRLB  is 
var(0)  =  cr^PDOP,  where  cr^  is  the  variance  of  the  range  measurement.  The  PDOP  is 
important  in  satellite  based  geolocation  especially  because  there  is  a  set,  predictable  orbit 
which  the  sensors  follow.  The  orbit  is  essentially  fixed  unless  the  satellites  are  equipped 
with  enough  fuel  to  perform  more  than  station  keeping  movement. 
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3.7  Simulation  Variables 


The  primary  drivers  of  the  results  of  the  simulator  described  in  this  chapter  are  the 
variables  which  modify  the  environment.  This  simulation  uses  four  variables  to  generate 
emitter  position  estimates.  First,  the  number  of  satellites  in  the  geolocation  constellation  is 
varied  from  one  to  four.  This  necessitates  the  use  of  different  algorithms  for  different 
configurations.  Some  algorithms,  such  as  TDOA  and  DPD,  cannot  work  with  a  single 
satellite.  Second,  the  SNR  of  the  received  signal  is  varied  from  -35  to  -l-35  dB  in  order  to 
simulate  different  transmitter  power  levels.  The  SNR  is  the  only  metric  used  in  this 
research  to  measure  the  quality  of  the  received  signal  itself.  Third,  the  size  of  the 
bounding  box  applied  to  the  target  search  grid  is  varied.  This  parameter  is  not  applicable 
in  some  algorithms,  such  as  AOA.  The  performance  of  the  algorithms  when  these 
variables  are  changed  provides  information  as  to  their  usefulness  in  different  scenarios, 
and  should  aid  designers  in  the  creation  of  effective  geolocation  systems. 

3.8  Summary 

This  chapter  has  provided  the  reader  with  a  detailed  description  of  the  simulation 
environment,  variables,  algorithms,  and  performance  metrics  used.  The  simulator 
presented  makes  use  of  a  robust  model  of  the  environment  in  which  the  geolocation 
system  must  operate,  and  will  provide  the  user  with  an  estimate  of  the  emitter’s  position 
for  a  given  set  of  environmental  variables. 

Chapter  IV  showcases  the  results  of  simulation,  and  indicates  how  the  chosen 
variables  will  affect  the  performance  of  the  geolocation  system. 
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IV.  Results  and  Analysis 


This  chapter  presents  analysis  of  the  results  of  the  simulations  described  in 

Chapter  III.  Multiple  scenarios  are  presented,  with  parameters  described  and  results 
analyzed  and  discussed.  The  performance  metrics  introduced  in  Chapter  II  are  analyzed  as 
well,  with  respect  to  the  scenarios  established  in  this  chapter.  Section  4. 1  introduces  the 
variables  common  to  all  of  the  scenarios  presented.  Sections  4.2  to  4.5  discuss  the 
scenarios  used  to  determine  geolocation  performance.  Sections  4.6  to  4.9  analyze  how 
some  variables  will  affect  the  system’s  geolocation  estimate.  Section  4.10  summarizes  the 
results  from  this  research. 

4.1  Global  Parameters 

The  global  parameters  used  in  simulation  are  those  which  are  true  for  the  entire 
ensemble  of  scenarios  presented  in  this  research.  The  SNR  used  in  each  scenario  is  varied 
from  -35  dB  to  -1-35  dB  in  10  dB  increments.  The  orbit  parameters  noted  in  Section  3.2  are 
not  changed  throughout  the  scenarios,  though  each  scenario  uses  a  different  number  of 
satellites.  A  total  of  515  discrete  1  ms  long  received  signals  are  used,  representing  the 
signal  collected  during  a  single  pass  over  the  target.  Timing  errors  associated  with  the 
satellites’  individual  ‘true’  time  result  from  use  of  the  GPS  system,  and  are  Gaussian  with 
zero  mean  and  a  standard  deviation  of  40  ns,  as  described  in  Section  3.3.  The  satellite’s 
knowledge  of  its  own  position  is  assumed  to  be  accurate  to  within  10  meters,  which  is 
typical  for  a  satellite  using  a  star  tracker  to  estimate  position.  The  simulated  emitter  is 
located  at  the  Air  Force  Institute  of  Technology,  and  is  stationary. 

4.2  Scenario  1  -  A  Single  Satellite 

A  single  satellite  relies  upon  a  changing  position  along  its  orbit  to  locate  a  stationary 
emitter  in  three  dimensions,  as  discussed  in  Section  2.1.1.  Simulations  for  the  single 
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satellite  seenario  use  the  AOA  and  IRF  methods.  These  two  methods  are  employed 
beeause  they  do  not  require  multiple  measurements  within  the  same  time  epoeh  to  solve 
for  the  emitter’s  position. 

The  position  estimate  is  extrapolated  from  lines  of  bearing  for  the  AOA  method  and  a 
grid  seareh  of  potential  emitter  loeations  for  the  IRF  method,  as  detailed  in  Seetions  3.5.2 
and  3.5.5,  respeetively.  The  IRF  method  measures  the  Doppler  of  the  reeeived  signal  over 
time,  and  eompares  the  measured  frequeney  to  the  ealeulated  frequeney  from  eaeh  grid 
point.  The  true  frequeney  and  simulated  instantaneous  reeeived  frequeney  at  an  SNR  of 
+5  dB  is  shown  in  Figure  4.1.  Notiee  the  ‘bunehing  up’  of  the  measurements  at  the  edges 
of  the  figure.  The  satellite’s  elevation  angle,  as  seen  from  the  transmitter,  will  ehange  more 
rapidly  as  the  satellite  reaehes  the  elosest  point  to  the  target  along  its  trajeetory  overhead. 


Doppler  Shift  Through  Single  Orbit 


Figure  4.1:  True  Doppler  and  Measured  Doppler  (SNR  =  +5  dB) 


In  this  seenario,  the  sensor  travels  along  an  orbit  aeeording  to  Table  3.1.  The 
simulation  is  eonsidered  with  20  trials  per  SNR  step.  Figure  4.2  shows  the  modeled 
performanee  of  the  AOA  and  IRF  methods  in  this  seenario. 
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Geolocation  Performance  for  One  Satellite 


Figure  4.2:  Algorithm  performance  for  a  single  satellite  scenario 


The  IRF  method  appears  to  perform  better  than  the  AOA  method  at  lower  SNR 
values,  and  the  AOA  method’s  performance  begins  to  diverge  from  the  IRF  method’s 
performance  when  the  SNR  reaches  -5  dB.  Both  estimators  continue  to  track  alongside  the 
CRLB  as  the  simulated  signal’s  SNR  increases.  However,  the  AOA  method  outperforms 
the  IRF  method  by  a  wide  margin  at  higher  SNR  values. 

The  trade  space  for  a  single  satellite  geolocation  system  design  includes  the  hardware 
requirements,  data  storage  and  throughput  requirements,  processing  requirements,  and 
performance  requirements.  A  significant  difference  in  the  hardware  requirements  of  these 
two  methods  is  the  number  of  antennas  required  to  perform  the  relevant  measurements. 
The  AOA  method  requires  a  minimum  of  four  antennas  to  achieve  a  three  dimensional 
bearing  estimate  (6,  (p).  More  antennas  improve  the  accuracy  of  the  MUSIC  estimation 
algorithm,  so  the  antenna  array  design  for  a  Cubesat  employing  the  AOA  method  is  an 
important  part  of  the  overall  system  design.  The  IRF  method  requires  a  single  antenna,  as 
the  only  measurement  is  of  the  instantaneous  frequency  of  the  received  signal.  The 
hardware  requirements  are  not  a  focus  of  this  research,  but  are  mentioned  to  give  the 


65 


reader  an  idea  of  the  trade  space.  The  IRF  method  is  much  more  sensitive  to  the  duration 
of  a  received  signal  than  the  AOA  method  is,  so  it  will  require  more  data  storage  and 
handling.  The  performance  requirements  depend  greatly  on  the  expected  signal  strength, 
hardware  available,  orbit  parameters,  and  other  factors.  If  the  best  possible  performance  is 
desired,  the  AOA  method  will  generally  perform  better  than  the  IRF  method.  However,  if 
hardware  space  is  limited  and  a  slightly  diminished  geolocation  performance  is 
acceptable,  the  IRF  method  may  be  a  better  option.  The  Cubesat  architecture  is  small 
enough  that  the  value  of  less  or  smaller  hardware  is  often  higher  than  the  performance  of 
the  payload,  so  the  IRF  method  may  be  a  better  option,  depending  on  the  circumstances. 

4.3  Scenario  2  -  Two  Satellites 

Adding  a  second  satellite  increases  the  performance  of  the  geolocation  system.  A 
second  satellite  provides  twice  the  information  per  time  epoch,  allowing  each  method  to 
converge  to  a  solution  faster.  Although  adding  a  second  satellite  inevitably  induces  a  time 
synchronization  error  between  the  two  satellites,  other  errors  are  reduced.  Two  satellites 
offer  a  dimensional  aspect  to  the  receivers  that  was  not  present  with  a  single  satellite.  The 
satellites  maintain  a  distance  of  approximately  70  km  from  each  other  while  in  the  range 
of  the  emitter’s  signal.  The  angular  separation  of  the  two  satellites  as  seen  from  the 
transmitter’s  point  of  view  is  at  most  9°.  While  this  separation  is  very  small,  the  slight 
offset  in  orbit  between  the  two  receivers  will  provide  small  differences  in  the  Doppler 
frequency  and  angle  measurements  made  by  the  platforms.  These  differences  help 
distinguish  the  closest  grid  point  in  the  IRF  method,  and  will  provide  the  AOA  method 
with  more  LOB’s  to  use  for  its  localization  process. 

In  this  scenario,  the  sensors  travel  along  orbits  as  described  in  Table  3.1.  The 
simulation  was  run  with  20  trials  per  SNR  step,  and  converges  toward  the  CRLB  as  the 
SNR  is  increased.  Figure  4.3  illustrates  the  simulated  performance  of  the  AOA  and  IRF 
methods  in  this  scenario. 
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Geolocation  Performance  for  Two  Satellites 


Figure  4.3:  Algorithm  performance  for  a  dual  satellite  scenario. 


The  AOA  method  performs  better  than  the  IRF  method  at  higher  SNR  values,  similar 
to  the  single  satellite  scenario.  However,  the  gap  in  performance  between  the  two  methods 
is  smaller  when  two  satellites  are  present. 

The  two  satellite  scenario  must  include  trade  space  considerations,  similarly  to  the 
single  satellite  scenario,  but  will  add  some  new  conditions.  If  two  satellites  are  to  work 
cooperatively,  a  communications  link  must  exist  between  satellites  if  any  joint  processing 
is  to  be  done  on-board.  The  additional  communications  hardware  will  take  up  space, 
which  contributes  to  the  argument  against  the  larger  payload  hardware  that  the  AOA 
method  requires.  The  data  that  must  be  transferred  between  the  two  sensor  platforms  is 
relatively  small  for  each  method.  The  AOA  method  only  requires  the  transmission  of  the 
fob’s  calculated  and  the  satellite  positions  for  each  time  epoch.  The  IRF  method  requires 
the  transmission  of  the  measured  instantaneous  received  frequency  and  the  satellite’s 
position  for  each  time  epoch.  Again,  the  size  of  the  platform  and  the  good  performance  of 
the  IRF  method  may  make  it  a  well  balanced  choice  for  system  designers. 
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4.4  Scenario  3  -  Three  Satellites 


A  third  satellite  enables  the  use  of  DPD  and  TDOA.  The  DPD  method  makes  use  of 
both  time  and  frequency  information,  as  discussed  in  Section  3.5.4.  TDOA  can  be 
performed  with  the  constraint  that  the  emitter’s  position  is  on  the  surface  of  the  Earth,  as 
detailed  in  Section  2.3. 3. 2.  The  AOA  and  IRE  methods  are  still  used  in  simulations  to 
compare  performance,  but  in  this  scenario  with  three  measurements  per  time  epoch.  A 
constellation  of  three  satellites  will  also  permit  the  geometry  of  the  receivers  in  a  triangle 
formation  rather  than  a  line,  which  is  the  case  when  using  just  two  satellites.  The  triangle 
geometry  improves  the  geolocation  estimate,  and  the  PDOP  may  now  be  calculated. 

The  geometry  between  three  satellites  and  the  transmitter  is  much  improved  over 
having  just  two  satellites.  The  PDOP  improves  as  the  satellites  move  across  the  target 
overhead.  Also,  the  orbit  parameters  for  the  third  satellite  are  adjusted  in  inclination  so 
that  the  constellation  is  at  its  widest  while  over  the  target.  The  PDOP  reaches  its  minimum 
of  53  when  the  satellites  are  directly  over  the  target. 

The  performance  of  each  of  the  four  algorithms  is  shown  in  Eigure  4.4.  The  four 
algorithms  are  run  in  20  trials  over  the  range  of  SNR  values. 

The  DPD  method  outperforms  the  AOA  and  IRE  methods  by  more  than  an  order  of 
magnitude  for  the  span  of  -15  to  +15  dB  SNR.  TDOA  nearly  matches  the  performance  of 
DPD  above  -5  dB  SNR,  and  all  methods  approach  the  CREB  as  the  SNR  is  increased.  The 
AOA  and  IRE  methods  both  improve  performance  with  three  satellites,  but  do  not  match 
the  DPD  or  TDOA  methods  in  performance.  It  is  interesting  to  note  that  the  TDOA 
method  outperforms  the  DPD  method  at  both  extremes  of  the  SNR  used  in  the  simulation. 
The  bounding  box  for  the  DPD  method  plays  a  part  in  the  decreased  performance  at  low 
SNRs.  The  bounding  box  for  this  simulation  is  2,000  km  per  side.  The  simulation  shows 
that  the  performance  of  the  DPD  method  at  less  than  -25  dB  is  highly  dependent  on  the 
size  of  the  emitter  location  grid  used  in  the  DPD  algorithm.  Eor  high  SNR,  the  TDOA 
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Geolocation  Performance  for  Three  Satellites 


Figure  4.4:  Performance  of  geolocation  methods  with  a  three-satellite  constellation 


method  again  outperforms  the  DPD  method.  The  TDOA  algorithm  for  three  satellites 
includes  the  constraint  that  the  transmitter  be  located  on  the  Earth’s  surface.  The  DPD 
method  does  not  require  such  a  constraint,  so  the  TDOA  method  ultimately  has  more 
information  than  the  DPD  algorithm,  leading  to  better  performance.  The  DPD  method 
consistently  outperforms  the  TDOA  method  for  a  large  range  of  SNR,  and  does  so  without 
constraints  on  the  altitude  of  the  emitter.  If  the  emitter  does  not  lie  on  the  surface  of  the 
earth,  then  the  TDOA  method  cannot  reliably  be  used  with  only  three  satellites. 

The  trade  space  for  a  three  satellite  system  has  the  same  considerations  as  for  two 
satellites.  The  algorithms  that  can  be  used  by  expanding  from  two  to  three  satellites  have 
implications  for  the  choice  of  methods  to  use  in  the  system.  The  DPD  and  TDOA  methods 
only  require  a  single  receiving  antenna  per  sensor,  with  an  additional  cross-link  antenna 
required  if  processing  is  completed  on-board.  The  computational  power  required  for  the 
DPD  method  is  similar  to  the  IRE  method;  a  grid  of  CAE  values  must  be  populated  and 
exhaustively  searched  for  each  possible  emitter  position  at  each  time  epoch.  The 
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computational  requirements  for  the  DPD  method  are  typically  greater  than  that  of  the 
TDOA,  where  a  simple  cross  correlation  between  signals  is  sufficient  to  determine  timing 
information.  Data  transfer  for  Weiss’  DPD  method  would  also  be  larger  than  other 
methods  because  a  representation  of  the  entire  signal  for  each  collection  must  be 
rebroadcast  to  a  single  computing  element.  TDOA  only  requires  the  transfer  of  timing 
information  between  sensors,  which  is  a  relatively  small  transfer.  Pourhomayoun  and 
Fowler  propose  methods  [16]  to  decrease  computational  complexity  and  data  transfer 
requirements  for  Weiss’  DPD  method.  This  research  did  not  implement  these  proposed 
improvements,  and  further  research  into  the  implications  of  such  modifications  for  a 
Cubesat-based  geolocation  mission  would  be  appropriate.  With  the  above  considerations 
in  mind,  the  DPD  method  has  the  best  performance  of  the  methods  outlined,  and  has  the 
potential  to  be  an  excellent  method  of  geolocation  for  a  space  platform. 

4.5  Scenario  4  -  Four  Satellites 

A  four  satellite  geolocation  system  allows  the  full  use  of  TDOA  without  constraints. 
Additionally,  it  creates  an  over-determined  solution  for  AOA,  IRF,  and  DPD,  which  can 
improve  the  solution.  The  four  satellites  are  configured  in  a  parallelogram  formation, 
which  changes  size  and  shape  as  the  orbit  progresses.  The  increase  in  sensors  will 
generally  decrease  the  PDOP  for  the  solution,  as  shown  in  Figure  4.5. 

The  PDOP  curve  for  the  four-satellite  constellation  is  similar  to  the  three-satellite 
constellation,  but  achieves  a  lower  minimum  of  41.  In  order  to  improve  the  PDOP  further, 
the  distance  between  the  satellites  could  be  increased,  improving  the  geometry  of  the 
constellation  over  the  transmitter. 

The  performance  of  the  algorithms  is  shown  in  Figure  4.6.  The  four  algorithms  are 
run  in  20  trials  over  the  range  of  SNR  values. 
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PDOP  for  4  Satellites 


Figure  4.5:  PDOP  for  a  four  satellite  eonstellation  in  orbit  over  the  target  emitter 


Geolocation  Performance  for  Four  Satellites 


Figure  4.6:  Performanee  of  geoloeation  methods  with  a  four-satellite  eonstellation 


The  performanee  of  the  DPD  method  with  four  satellites  in  the  eonstellation  is  better 
than  the  other  methods  by  a  large  margin.  The  SNR  range  between  -35  dB  and  -1-15  dB 
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again  shows  that  the  DPD  method  excels  compared  to  other  methods  when  the  signal 
strength  is  low.  The  size  of  the  bounding  box  used  in  the  grid  search  methods  was  200  km 
per  side,  which  could  model  a  small  country.  The  performance  of  the  TDOA  method 
actually  decreases  in  this  scenario  because  the  constraint  that  the  emitter  be  on  the  Earth’s 
surface  is  removed.  Without  that  constraint,  geolocation  accuracy  using  the  TDOA 
method  is  diminished  because  the  range  to  the  transmitter  is  so  much  larger  than  the  range 
between  satellites.  If  the  mission  designer  knew  the  transmitter  was  on  the  surface  of  the 
Earth,  the  constraint  could  remain  and  the  performance  would  likely  improve  slightly 
from  the  three  satellite  scenario. 

The  trade  space  considerations  for  the  four  satellite  system  are  unchanged  from  the 
three  satellite  system.  Again,  DPD  appears  to  offer  the  best  performance,  and  may  have 
the  best  value  to  a  Cubesat-based  geolocation  system  with  a  constellation  of  three  or  more 
satellites.  The  drawback  of  this  method  is  the  requirement  that  the  signal  is  transmitted  to 
a  central  location  for  processing.  An  example  signal  composed  of  5000  data  points  in 
Matlab  creates  a  file  which  is  30  kilobytes.  A  collection  of  515  signal  samples  taken 
during  as  pass  over  the  target  would  be  approximately  15  megabytes  in  size. 

4.6  Effects  of  Increasing  Constellation  Size  for  the  AOA  Method 

The  effect  of  adding  more  satellites  to  a  geolocation  system  is  somewhat  obscured  by 
the  performance  plots  provided  in  Sections  4.2  to  4.5.  In  order  to  see  the  effect  of  adding 
satellites  to  the  constellation  more  clearly,  the  results  can  be  shown  as  a  function  of 
measurements  per  estimate  rather  than  SNR.  The  SNR  is  held  constant  at  +5  dB,  while  the 
number  of  satellites  is  varied  from  one  to  four.  The  AOA  algorithm  is  used  because  it  can 
be  emlpoyed  in  all  four  scenarios.  Eigure  4.7  shows  the  convergence  of  the  solution  given 
for  different  satellite  configurations. 
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Figure  4.7:  AOA  Error  shown  as  a  function  of  the  number  of  satellites  in  the  system 


The  error  shown  in  Figure  4.7  approaches  zero  as  more  measurements  are  obtained. 
Increasing  the  number  of  satellites  decreases  the  time  required  for  convergence,  and 
results  in  a  smaller  error  magnitude  when  convergence  is  achieved.  When  a  single  satellite 
is  present,  the  solution  doesn’t  really  converge  until  the  orbit  is  nearly  complete  at  515 
samples.  When  four  satellites  are  present,  the  solution  converges  before  the  orbit  is 
halfway  complete,  at  around  200  measurements.  This  is  because  there  is  four  times  as 
much  information  per  measurement  when  four  satellites  are  present,  compared  to  a  single 
satellite. 

4.7  Effects  of  Changing  the  Bounding  Box  on  DPD 

The  DPD  and  IRE  methods  use  a  grid  search  to  perform  the  position  calculation, 
using  a  predetermined  ‘bounding  box’  within  which  the  emitter  is  assumed  to  be  located. 
The  bounding  box  is  used  to  define  initial  grid  parameters,  and  the  box  becomes  smaller 
as  the  algorithm  progresses,  until  the  final  bounding  box  is  a  cube  with  10  meter  sides. 

The  initial  bounding  box  size  can  be  altered  in  the  simulation  to  show  its  effect  on  the 
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performance  of  the  locator.  Figure  4.8  shows  the  performance  of  the  DPD  method  using 
bounding  box  sizes  varying  from  100  m  to  1000  km.  Three  satellites  are  used  in  each  trial, 
with  10  trials  per  SNR  step.  The  same  SNR  step  of  10  dB  is  used  for  SNR  values  ranging 
from  -35  dB  to  -1-35  dB.  The  grids  are  designed  as  described  in  Section  3.5.4,  with  125 
points  in  the  grid,  and  decreasing  in  size  by  20  percent  per  dimension  in  each  iteration. 
The  results  are  shown  in  Figure  4.8. 


Figure  4.8:  The  performance  of  the  DPD  method  is  dependent  on  the  size  of  the  bounding 
box,  especially  at  lower  SNR  values 


Figure  4.8  shows  that  the  performance  of  the  DPD  method  is  dependent  on  the 
bounding  box  used  in  the  algorithm,  as  well  as  signal  strength.  The  results  show  that  the 
algorithm  provides  a  solution  at  the  edges  of  the  bounding  box  when  the  SNR  is  at  or 
below  -25  dB.  As  the  SNR  improves,  the  solution  will  improve  as  well.  The  figure  shows 
that  bounding  boxes  with  sides  smaller  than  10  km  have  a  bias  even  when  SNR  is  greater 
than  -15  dB.  The  bounding  box  must  be  taken  into  consideration  when  designing  a 
geolocation  mission  using  the  DPD  or  IRF  method  because  the  algorithm  relies  so  heavily 
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on  the  correct  placement  and  size  of  the  search  grid.  If  the  bounding  box  is  too  small  and 
is  misplaced,  the  method  may  converge  to  an  incorrect  solution  simply  because  the  true 
emitter  position  is  not  present  within  the  perimeter  of  the  search  grid.  Additionally,  the 
mission  designer  must  assess  whether  Cubesat  geolocation  provides  information  if  the 
bounding  box  is  small.  If  the  target  position  is  known  to  be  within  a  I  km  block,  the  DPD 
method  will  not  provide  better  information  than  this  unless  the  signal  SNR  at  the  receiver 
is  greater  than  +15  dB. 

4.8  CRLB  Analysis 

The  CRLB  for  each  parameter  is  calculated.  The  CRLB  for  bearing  measurements 
offers  interesting  results  for  the  case  of  a  fixed  orbit.  A  normal  orbit  has  a  period  of  time 
where  the  satellite  is  approaching  the  target  and  a  period  of  time  when  the  satellite  is 
moving  away  from  the  target.  The  lower  bound  of  the  bearing  estimation  accuracy  is 
highly  dependent  on  the  true  angle  from  which  the  signal  is  arriving.  The  bearing 
estimation  accuracy  is  maximized  when  the  signal  arrives  at  the  antenna  array  when  the 
satellite  reaches  zenith  relative  to  the  target  emitter.  The  bearing  estimate  accuracy 
reaches  a  minimum  when  the  signal  reaches  the  the  array  at  0°.  The  CRLB  for  the 
elevation  angle  illustrated  in  Figure  4.9  showcases  this  phenomenon  for  an  SNR  of  +15 
dB.  When  the  signal  is  initially  received  by  the  sensor  as  it  orbits,  the  target  is  presumably 
on  or  near  the  horizon.  When  the  satellite  reaches  the  smallest  distance  from  the  target, 
the  estimate  will  be  much  more  reliable  because  the  true  elevation  angle  will  be  higher. 
Figure  3.3  showed  that  the  satellites  in  this  research  will  not  pass  directly  overhead  the 
target,  so  the  elevation  angle  will  not  reach  90°.  However,  Figure  4.9  shows  that  as  the 
satellite  approaches  the  target,  the  CRLB  will  decrease,  and  as  the  satellite  moves  away 
from  the  target,  the  CRLB  will  increase  again. 
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Cramer-Rao  Lower  Bound:  Bearing  -  Elevation 


Figure  4.9:  The  elevation  angle  CRLB  calculated  over  one  pass  will  decrease  as  the  true 
bearing  angle  increases  toward  90°,  then  increase  again  as  the  bearing  angle  decreases 
toward  0° 


This  research  assumes  that  one  set  of  sensing  antennas  is  placed  facing  nadir,  and  one 
set  is  placed  90°  offset  from  nadir.  Additionally,  the  satellite  is  assumed  to  have  enough 
fuel  to  perform  station-keeping  maneuvers  throughout  its  lifetime.  The  azimuth  bearing 
for  this  particular  scenario  never  reaches  0°,  and  fares  much  better  than  the  elevation  angle 
throughout  the  orbit,  as  shown  in  Figure  4.10.  The  orientation  of  the  antennas  on  the 
satellite’s  frame  will  have  an  effect  on  the  angle  at  which  the  signal  is  sensed  as  well,  and 
could  change  this  CRLB  significantly. 

The  frequency  and  timing  lower  bounds  do  not  change  with  target-to-satellite 
geometry  as  with  the  lower  bound  on  bearing  angles.  As  detailed  in  Section  3. 6. 1.1,  the 
lower  bounds  for  frequency  and  timing  are  dependent  on  signal  parameters  such  as  SNR 
and  bandwidth.  Figure  4.1 1  and  Figure  4.12  show  the  CRLB  on  standard  deviation 
calculated  for  the  frequency  measurement  and  timing  measurement,  respectively. 
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Cramer-Rao  Lower  Bound:  Bearing  -  Azimuth 


Figure  4.10:  CRLB  calculated  for  azimuth  angle  estimation  over  one  orbit 


Figure  4.1 1:  CRLB  calculated  for  frequency  measurement  at  a  range  of  SNR  values 


Figure  4.11  shows  that  the  frequency  measurement  for  this  simulation  has  a  statistical 
lower  bound  of  about  200  Hz  at  an  SNR  of  +35  dB.  The  timing  measurement  has  a 
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Figure  4.12:  CRLB  calculated  for  timing  measurement  at  a  range  of  SNR  values 


statistical  lower  bound  of  about  40  microseconds  at  an  SNR  of  +35  dB,  as  in  Figure  4.12. 
These  figures  could  be  improved  by  increasing  the  bandwidth  of  the  signal  and  the 
duration  of  the  collection. 

4.9  PDOP  Results 

The  PDOP  calculated  for  the  satellite  constellation  configurations  used  in  this 
research  shows  the  effect  geometry  has  on  the  performance  of  the  algorithms.  As  the 
results  in  the  previous  section  showed,  adding  more  satellites  will  promote  a  faster 
convergence  to  a  geolocation  solution.  Figure  4.13  shows  the  PDOP  calculated  for  the 
scenarios  with  three  and  four  satellite  constellations  (PDOP  cannot  be  reasonably 
computed  for  fewer  than  three  satellites). 

The  PDOP  is  shown  to  be  a  curve,  where  the  orbit  has  the  lowest  PDOP  at  zenith. 
The  best  geometry  is  at  the  center  of  the  orbit  because  the  ratio  between  the  target  range 
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PDOP  for  3  and  4  Satellite  Constellations 


Figure  4.13:  Adding  satellites  to  the  geoloeation  eonstellation  is  shown  to  improve  PDOP 


and  the  satellite  separation  is  at  its  smallest.  According  to  Langley,  DOP  will  always 
decrease  when  more  sensors  are  present  [42]. 

By  way  of  comparison,  the  PDOP  is  poor  for  all  of  the  scenarios  discussed  in  this 
chapter  compared  to  the  PDOP  that  a  typical  GPS  receiver  has.  GPS  typically  has  DOP 
measurements  less  than  10,  whereas  this  system  has  a  minimum  PDOP  of  41.  There  are 
elements  of  the  design  that  favor  smaller  distances  between  satellites,  which  will  increase 
the  PDOP.  For  example,  communication  links  between  satellites  are  more  effective  at 
smaller  distances,  and  could  constrain  the  geometry.  A  smaller  constellation  geometry 
may  also  increase  the  time  in  which  all  the  satellites  receive  the  SOI.  The  orbital  dynamics 
of  the  geoloeation  system  must  be  carefully  planned  to  balance  the  needs  of  the  satellite 
constellation  with  the  desire  to  minimize  the  dilution  of  precision  and  thereby  improve  the 
position  estimate. 
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4.10  Summary 

The  results  and  analysis  presented  in  this  chapter  characterizes  the  performance  of 
four  geolocation  algorithms  applied  to  a  number  of  sensors  orbiting  over  an  RF  emitter. 
The  analysis  of  the  scenarios  shows  that  the  DPD  method  generally  performs  better  than 
the  IRF  method,  as  well  as  the  classical  methods  of  TDOA  and  AOA.  The  performance 
gained  by  using  DPD  at  lower  SNR  values,  especially  in  the  region  from  -15  to  +15  dB,  is 
significant.  A  geolocation  mission  using  Cubesats  would  gain  a  substantial  improvement 
in  performance  by  implementing  the  DPD  algorithm. 

Several  considerations  for  the  selection  of  a  geolocation  method  are  addressed  in  this 
research.  The  hardware  requirements  for  the  AOA  method  include  an  antenna  array  rather 
than  a  single  antenna,  which  is  a  burden  for  such  a  small  satellite  architecture.  When 
multiple  satellites  are  present,  a  cross-link  radio  and  antenna  are  required  to  share 
information  between  the  satellites.  The  DPD  method  presented  by  Weiss  is  more 
computationally  complex  and  requires  more  data  transfer  than  the  other  methods 
considered,  but  these  issues  may  be  diminished  by  implementing  the  approximated  DPD 
method  proposed  by  Pourhomayoun  and  Fowler  [16].  This  research  did  not  use  the 
approximated  DPD  method,  and  further  investigation  into  its  performance  from  an 
orbiting  sensor  platform  should  be  conducted.  Given  the  results  presented  in  Sections  4.2 
to  4.5,  the  IRF  method  is  suggested  for  consideration  for  single  satellite  and  two  satellite 
constellations.  The  DPD  method  is  recommended  for  use,  versus  TDOA,  for 
constellations  with  three  or  more  satellites. 

The  results  presented  in  this  thesis  are  comparable  with  other  results  in  open 
literature  [18],  [16].  The  algorithms  considered  in  this  research  are  shown  generally  to 
approach  the  CRLB,  the  statistical  lower  bound  on  the  performance  of  an  estimator.  Other 
research  has  shown  this  tendency,  but  this  research  is  specific  to  the  consideration  of 
sensors  in  orbit.  Standoff  distances  on  the  order  of  one  thousand  kilometers  have  not  been 
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previously  demonstrated  to  work  for  the  DPD  method.  This  thesis  is  believed  to  present 
the  first  such  results. 

Chapter  V  provides  conclusions  drawn  from  this  research  and  recommendations  for 
future  research. 
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V.  Conclusion 


5.1  Conclusions 

This  research  presents  a  simulation  environment  and  methodology  to  characterize  the 
geolocation  performance  of  a  constellation  of  nano-satellites  in  LEO.  The  simulation 
incorporated  realistic  orbits  created  in  AGI’s  STK  software  to  use  in  conjunction  with  four 
geolocation  algorithms.  The  results  of  the  simulation  and  the  analysis  presented  in  this 
research  provide  a  thorough  comparison  of  the  performance  that  should  be  expected  from 
the  algorithms  and  methods  used  in  a  variety  of  scenarios.  Other  considerations,  such  as 
statistical  lower  bounds,  the  effect  of  geometry,  and  the  number  of  satellites  in  a 
geolocation  system  are  addressed  with  respect  to  the  geolocation  methods  investigated  in 
this  research. 

The  evidence  supports  further  research  into  the  Direct  Position  Determination 
method,  as  it  appears  to  perform  better  than  other  methods  at  the  standoff  range  of  a 
satellite  system  with  three  or  more  platforms.  When  a  Cubesat-based  geolocation  system 
has  only  one  or  two  satellites,  the  IRE  method  provides  a  balance  between  performance 
and  hardware  requirements  that  make  it  attractive  for  use  with  the  nano-satellite 
architecture  in  some  geolocation  scenarios. 

A  Cubesat  based  geolocation  system  can  be  a  powerful  tool  to  locate  terrestrial 
communication  systems,  jammers,  radars,  and  other  RE  emitters.  This  research  has 
displayed  the  use  of  four  geolocation  methods  which  can  be  used  in  different 
circumstances  to  locate  an  emitter. 

When  a  single  Cubesat  is  used  to  geolocate,  the  AOA  and  IRE  algorithms  may  be 
used  because  the  satellite  has  velocity  while  moving  along  its  orbit.  The  AOA  and  IRE 
methods’  performance  suggest  that  the  AOA  method  will  generally  perform  better  on 
orbit.  While  the  AOA  method  outperforms  the  IRE  method,  the  smaller  footprint  of  the 
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hardware  requirements  for  the  IRF  method  make  it  appealing  if  some  performance 
degradation  is  acceptable.  By  adding  a  second  satellite,  the  performance  of  the  AOA  and 
IRF  methods  will  improve,  decreasing  the  time  required  to  approach  the  CRLB.  The 
evidence  supports  further  research  into  the  Direct  Position  Determination  method,  as  it 
appears  to  perform  better  than  other  methods  at  the  standoff  range  of  a  satellite  system 
with  three  or  more  platforms.  The  DPD  method  has  been  shown  to  be  effective  in 
geolocating  RF  emitters  at  satellite-range  standoff  distances.  When  a  Cubesat-based 
geolocation  system  has  only  one  or  two  satellites,  the  IRF  method  provides  a  balance 
between  performance  and  hardware  requirements  that  make  it  attractive  for  use  with  the 
nano-satellite  architecture. 

5.2  Recommendations  for  Future  Work 

There  are  a  number  of  areas  in  which  this  research  should  be  continued  and 
enhanced.  Including  a  robust  model  for  the  attenuation  associated  with  signal  broadcast 
through  the  atmosphere  would  make  the  propagation  of  the  transmitted  signal  more 
realistic.  The  atmospheric  attenuation  a  signal  experiences  is  dependent  on  the  frequency 
of  that  signal,  so  an  in-depth  analysis  of  the  impact  on  geolocation  performance  associated 
with  variation  in  the  SOI  carrier  frequency  should  be  completed  with  use  of  this  more 
robust  model. 

A  study  of  the  orbital  properties  of  the  system,  and  how  they  would  affect  the 
performance  of  geolocation  would  be  appropriate  as  well.  PDOP  would  be  an  important 
consideration  in  such  an  analysis.  Additionally,  more  passes  could  be  added  to  the 
simulation  to  show  how  more  collected  data  will  improve  the  solution. 

Further  improvements  can  be  made  to  the  algorithms  used  in  the  geolocation  system. 
The  IRF  method  would  benefit  from  improving  the  precision  and  processing  time  for 
measuring  the  instantaneous  received  frequency.  The  approximated  DPD  method  [16] 
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should  be  explored  in  detail  and  compared  to  the  performance  of  Weiss’  DPD  method  [18] 
to  ensure  performance  degradation  at  standoff  distances  is  sufficiently  small. 

The  type  of  signal  could  be  analyzed  with  respect  to  geolocation  performance.  This 
research  focused  on  a  narrowband  signal  with  a  relatively  generic  pulse  for  the  signal’s 
structure.  Wideband  and  ultra-wideband  signal  sources  are  common  in  the 
communications  field,  and  could  present  unique  challenges  for  geolocation  systems. 
Additionally,  different  modulation  techniques  such  as  Gaussian  minimum  shift  keying 
(OMSK)  or  binary  phase  shift  keying  (BPSK)  may  be  important  factors  in  the  design  of  a 
geolocation  system. 
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